{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "fffbbaec",
   "metadata": {},
   "source": [
    "# 利用 Product Formula 模拟时间演化\n",
    "<em> Copyright (c) 2021 Institute for Quantum Computing, Baidu Inc. All Rights Reserved. </em>"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e6f9a147",
   "metadata": {},
   "source": [
    "## 概述\n",
    "\n",
    "量子力学中系统的能量由哈密顿量算符 $H$ 描述，求解给定系统的哈密顿量的全部或者部分性质，构成了凝聚态物理、计算化学和高能物理等一系列学科的核心问题。然而，由于系统的自由度随系统增大呈指数级增加，导致一般情况下无法利用经典计算机有效模拟量子系统——即便用光全世界的内存和硬盘也不能直接储存一个仅几百量子比特的系统的态矢量。与经典计算机不同，量子计算机由于其所有的操作都是直接作用在同样为指数级别的态空间上，因而在模拟一个量子系统上量子计算机具有经典计算机无法比拟的优势。实际上，设计一个可控的量子系统来高效模拟自然界中的量子系统，正是费曼在上世纪 80 年代提出量子计算这一概念时的最初想法：\n",
    " \n",
    "> _\"Nature isn't classical, dammit, and if you want to make a simulation of nature, you'd better make it quantum mechanical, and by golly it's a wonderful problem, because it doesn't look so easy.\"_\n",
    ">\n",
    "> --- \"Simulating physics with computers\", 1982, Richard P. Feynman [1]\n",
    "\n",
    "通用量子计算机以及一系列量子模拟器的发展令费曼的设想有了实现的可能。在通用量子计算机上进行数字量子模拟（digital quantum simulation）—— 利用量子门构造量子线路从而实现量子模拟，由于该方法具有较高的可拓展性和通用性，因而被认为是最有潜力的技术路线。\n",
    "\n",
    "本教程讲述了如何利用 Paddle Quantum 模拟量子系统的时间演化过程，主要分为以下三个部分：\n",
    "1. 如何在 Paddle Quantum 中创建并操作一个 `Hamiltonian` 对象\n",
    "2. 如何利用 `construct_trotter_circuit` 来构建时间演化电路\n",
    "3. Suzuki product formula 方法的原理，以及如何进一步搭建任意阶的 Trotter-Suzuki 电路\n",
    "\n",
    "\n",
    "## 定义系统哈密顿量\n",
    "\n",
    "在介绍如何搭建时间演化电路之前，读者可以先熟悉一下 Paddle Quantum 中的哈密顿量。目前用户需要通过定义一个列表的形式来创建哈密顿量，列表中的元素应为哈密顿量中的每一项的系数与其泡利算符。作为教程，我们首先考虑一个简单的哈密顿量：\n",
    "\n",
    "$$\n",
    "H = Z \\otimes Z\n",
    "\\tag{1}\n",
    "$$\n",
    "\n",
    "该哈密顿量描述了两个量子比特之间的一种简单相互作用：当两个量子比特同时处于 $|0\\rangle$ 态或者 $|1\\rangle$ 态时，系统的能量为 $+1$；相反，当两个量子比特的态不同时，系统的能量为 $-1$。\n",
    "\n",
    "接下来，用户可以在 Paddle Quantum 中创建这一哈密顿量："
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "id": "1a17d7d8",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "1.0 Z0, Z1\n"
     ]
    }
   ],
   "source": [
    "import warnings\n",
    "warnings.filterwarnings(\"ignore\")\n",
    "from paddle_quantum.hamiltonian import Hamiltonian\n",
    "\n",
    "h = Hamiltonian([[1, 'Z0, Z1']])\n",
    "print(h)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "64ef6a63",
   "metadata": {},
   "source": [
    "目前，Paddle Quantum 中的哈密顿量类 `Hamiltonian` 支持自动合并同类项，加减法、索引以及拆分等操作："
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "ff08a2a7",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "1.0 Z0, Z1\n"
     ]
    }
   ],
   "source": [
    "h = Hamiltonian([[0.5, 'Z0, Z1'], [0.5, 'Z1, Z0']], compress=True)\n",
    "print(h)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "fec891a5",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "h + h: 2.0 Z0, Z1\n",
      "h * 2: 2.0 Z0, Z1\n",
      "h: 1.0 Z0, Z1\n"
     ]
    }
   ],
   "source": [
    "print('h + h:', h + h)\n",
    "print('h * 2:', h * 2)\n",
    "print('h:', h[:])"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d2be41b1",
   "metadata": {},
   "source": [
    "同时，内置的 `decompose_pauli_words()` 和 `decompose_with_sites()` 方法可以将哈密顿量分解为更加方便处理的形式："
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "40fcc0b6",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Pauli words 分解： ([1.0], ['ZZ'])\n",
      "Pauli with sites 分解： ([1.0], ['ZZ'], [[0, 1]])\n"
     ]
    }
   ],
   "source": [
    "print('Pauli words 分解：', h.decompose_pauli_words())\n",
    "print('Pauli with sites 分解：', h.decompose_with_sites())"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "95945ddd",
   "metadata": {},
   "source": [
    "除此之外，`construct_h_matrix()` 还可以创建其在泡利 $Z$ 基底下的矩阵形式："
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "497ff8fd",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "matrix([[ 1.+0.j,  0.+0.j,  0.+0.j,  0.+0.j],\n",
       "        [ 0.+0.j, -1.+0.j,  0.+0.j,  0.+0.j],\n",
       "        [ 0.+0.j,  0.+0.j, -1.+0.j,  0.+0.j],\n",
       "        [ 0.+0.j,  0.+0.j,  0.+0.j,  1.+0.j]], dtype=complex64)"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "h.construct_h_matrix()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "719696dd",
   "metadata": {},
   "source": [
    "## 模拟时间演化\n",
    "\n",
    "根据量子力学的基本公理，在确定了一个系统的哈密顿量之后，该系统随时间演化的过程可以由如下方程描述\n",
    "\n",
    "$$\n",
    "i \\hbar \\frac{\\partial}{\\partial t} | \\psi \\rangle = H | \\psi \\rangle,\n",
    "\\tag{2}\n",
    "$$\n",
    "\n",
    "$\\hbar$ 为约化普朗克常数。该方程正是著名的薛定谔方程。因此，对于一个不含时的哈密顿量，系统的时间演化方程可以写为\n",
    "\n",
    "$$\n",
    "|\\psi(t) \\rangle = U(t) | \\psi (0) \\rangle, ~ U(t) = e^{- i H t}.\n",
    "\\tag{3}\n",
    "$$\n",
    "\n",
    "这里我们取自然单位 $\\hbar=1$，$U(t)$ 为时间演化算符。利用量子线路来模拟时间演化过程的核心思想是利用量子电路构建出的酉变换模拟和近似该时间演化算符。在量桨中，我们提供了 `construct_trotter_circuit(circuit, Hamiltonian)` 函数来构建对应某一个哈密顿量的时间演化电路。下面，就让我们用刚刚创建的哈密顿量来测试一下："
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "3449029b",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "--*-----------------*--\n",
      "  |                 |  \n",
      "--x----Rz(2.000)----x--\n",
      "                       \n"
     ]
    }
   ],
   "source": [
    "from paddle_quantum.trotter import construct_trotter_circuit\n",
    "from paddle_quantum.ansatz import Circuit\n",
    "\n",
    "cir = Circuit()\n",
    "construct_trotter_circuit(cir, h, tau=1, steps=1) \n",
    "print(cir)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e3dc6c7f",
   "metadata": {},
   "source": [
    "可以看到，我们已经创建了针对 `h` 的基本量子模拟电路，它可以根据输入的 `tau` 来模拟任意时间长度的时间演化。\n",
    "\n",
    "此时，如果检查该电路对应的酉变换的矩阵形式，会发现它与时间演化算符 $e^{-iHt}$ 完全一致。这里，我们先使用 `gate_fidelity` 来计算量子电路的酉矩阵和时间演化算符的酉矩阵之间的保真度。保真度越接近 1 时代表两个酉矩阵代表的演化过程越相似。在下文，我们还会引入更加严格的误差定义。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "8d312985",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "电路的酉算符和时间演化算符之间的保真度为： 1.000\n"
     ]
    }
   ],
   "source": [
    "from scipy import linalg\n",
    "from paddle_quantum.qinfo import gate_fidelity\n",
    "\n",
    "# 计算 e^{-iHt} 和电路的酉矩阵之间的保真度\n",
    "print('电路的酉算符和时间演化算符之间的保真度为： %.3f' \n",
    "      % gate_fidelity(cir.unitary_matrix(), linalg.expm(-1 * 1j * h.construct_h_matrix())))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "fa34ba8a",
   "metadata": {},
   "source": [
    "实际上，这是因为泡利算符组成的张量积对应的任意角度旋转算符都可以被很好的分解为门电路的形式。比如在这个例子中我们只需要改变电路中 Rz 门的角度，就可以实现任意 $e^{-i Z\\otimes Z t}$ 演化算符的模拟。那么，这是否意味着可以用类似的电路形式来对任意的写成泡利形式的哈密顿量进行完美的模拟呢？答案是否定的。考虑一个稍微复杂一些的哈密顿量：\n",
    "\n",
    "$$\n",
    "H = Z \\otimes Z + X \\otimes I + I \\otimes X.\n",
    "\\tag{4}\n",
    "$$\n",
    "\n",
    "同样地，让我们用 `construct_trotter_circuit` 来构建它所对应的时间演化电路："
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "40bcf1e4",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "--*-----------------*----Rx(2.000)--\n",
      "  |                 |               \n",
      "--x----Rz(2.000)----x----Rx(2.000)--\n",
      "                                    \n",
      "电路的酉算符和时间演化算符之间的保真度为： 0.681\n"
     ]
    }
   ],
   "source": [
    "h_2 = Hamiltonian([[1, 'Z0, Z1'], [1, 'X0'], [1, 'X1']]) # 不需要写出单位算符\n",
    "cir = Circuit()\n",
    "construct_trotter_circuit(cir, h_2, tau=1, steps=1)\n",
    "print(cir)\n",
    "print('电路的酉算符和时间演化算符之间的保真度为： %.3f' \n",
    "      % gate_fidelity(cir.unitary_matrix().numpy(), linalg.expm(-1 * 1j * h_2.construct_h_matrix())))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "27f0a4ed",
   "metadata": {},
   "source": [
    "可以看到，此时电路对应的酉矩阵和时间演化矩阵之间的保真度 $<1$，说明该电路无法正确地模拟系统的时间演化过程。\n",
    "\n",
    "其实，对于我们此时构建的电路而言，它对应的酉变换为 $e^{-iZ\\otimes Z t} e^{-i (X\\otimes I + I\\otimes X)t}$，而时间演化算符为 $e^{-iZ\\otimes Z t - i(X\\otimes I + I\\otimes X)t}$。而对于一个量子系统而言，当两个算符不对易，即 $[A, B] \\neq 0$ 时，$e^{A+B} \\neq e^A e^B$。在这里，因为泡利算符 $X$ 和 $Z$ 之间并不对易，所以电路对应的酉变换并不等于正确的时间演化算符。\n",
    "\n",
    "除了利用保真度来描述量子电路和希望模拟的时间演化算符之间的相似性之外，我们还可以定义如下的误差 $\\epsilon$ 来刻画模拟演化电路和演化算符之间的距离，\n",
    "\n",
    "$$\n",
    "\\epsilon(U) = \\Vert e^{-iHt} - U \\Vert,\n",
    "\\tag{5}\n",
    "$$\n",
    "\n",
    "其中 $\\Vert \\cdot \\Vert$ 表示取最大的本征（奇异）值的模。这样的定义比起保真度而言更好的描述了量子态在不同的演化算符作用下产生的偏差，是一种更加严谨的定义模拟时间演化误差的方式。我们在下文中也会多次用到该误差。"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9a0f9969",
   "metadata": {},
   "source": [
    "### Product formula 与 Suzuki 分解\n",
    "\n",
    "Seth Lloyd 在 1996 年的文章中指出，可以将一整段的演化时间 $t$ 拆分为 $r$ 份较短的“时间块”来减小模拟时间演化的误差 [2]。考虑一个更一般的哈密顿量形式 $H = \\sum_{k=1}^{L} h_k$，其中 $h_k$ 是作用在一部分系统上的子哈密顿量。通过泰勒展开，不难发现其模拟误差是一个二阶项，即\n",
    "\n",
    "$$\n",
    "e^{-iHt} = \\prod_{k=1}^{L} e^{-i h_k t} + O(t^2).\n",
    "\\tag{6}\n",
    "$$\n",
    "\n",
    "那么，我们令 $\\tau = t/r$，并考虑演化算符 $\\left(e^{-iH \\tau}\\right)^r$，其演化误差为\n",
    "\n",
    "$$\n",
    "e^{-iHt} = \\left(e^{-iH \\tau}\\right)^r = \\left(\\prod_{k=1}^{L} e^{-i h_k \\tau} + O(\\tau^2) \\right)^r = \\left(\\prod_{k=1}^{L} e^{-i h_k \\tau} \\right)^r + O\\left(\\frac{t^2}{r}\\right).\n",
    "\\tag{7}\n",
    "$$\n",
    "\n",
    "上式告诉我们，只要将一整段演化时间拆为足够多的“片段”，就可以达到任意高的模拟精度，这就是 product formula 的基本思想。不过，(7) 中给出的误差只是一个粗略的估计。在实际情况中，为了估计达到某一模拟精度所需要的量子电路深度，就需要进一步推导其更严格的误差上界。下面，我们展示一个比较简略的误差上界计算过程，对具体的计算过程不感兴趣的读者可以直接跳至 (11) 阅读这部分的结论。\n",
    "\n",
    "记函数 $f$ 泰勒展开至 $k$ 阶之后的余项为 $\\mathcal{R}_k(f)$, 我们需要用到如下两个结论：\n",
    "\n",
    "$$\n",
    "\\left\\Vert \\mathcal{R}_k \\left( \\prod_{k=1}^{L} \\exp (-i h_k \\tau) \\right) \\right\\Vert\n",
    "\\leq\n",
    "\\mathcal{R}_k \\left( \\exp \\left( \\sum_{k=1}^{L} \\vert \\tau \\vert \\cdot \\Vert h_k \\Vert \\right) \\right),\n",
    "\\tag{8}\n",
    "$$\n",
    "\n",
    "$$\n",
    "\\vert \\mathcal{R}_k(\\exp (\\alpha)) \\vert \\leq \\frac{\\vert \\alpha \\vert^{k+1}}{(k+1)!} \\exp ( \\vert \\alpha \\vert ), ~\n",
    "\\forall \\alpha \\in \\mathbb{C}.\n",
    "\\tag{9}\n",
    "$$\n",
    "\n",
    "这两个结论的证明略微有一些繁琐，故在此略去，感兴趣的读者可以参考 [3] 中的 F.1 节。利用 (5) 中的定义，模拟误差可以写为\n",
    "\n",
    "$$\n",
    "\\epsilon\\left(e^{-iH\\tau}, U_{\\rm circuit}\\right) = \\left \\Vert \\exp\\left(-i\\sum_{k=1}^L h_k \\tau\\right) - \\prod_{k=1}^{L} \\exp\\left(-i h_k \\tau \\right) \\right \\Vert.\n",
    "\\tag{10}\n",
    "$$\n",
    "\n",
    "我们已经知道，此时的误差是泰勒展开至一阶后的余项，再利用 (8)，(9) 和三角不等式，可以将 (10) 中的误差上界进行计算：\n",
    "\n",
    "$$\n",
    "\\begin{aligned}\n",
    "\\left \\Vert \\exp\\left(-i\\sum_{k=1}^L h_k \\tau\\right) - \\prod_{k=1}^{L} \\exp\\left(-i h_k \\tau \\right) \\right \\Vert\n",
    "=~&\n",
    "\\left \\Vert \\mathcal{R}_1 \\left( \\exp \\left( -i \\sum_{k=1}^{L} h_k \\tau \\right) - \\prod_{k=1}^{L} \\exp (-i h_k \\tau) \\right) \\right \\Vert\n",
    "\\\\\n",
    "\\leq~&\n",
    "\\left \\Vert \\mathcal{R}_1 \\left( \\exp \\left( -i \\sum_{k=1}^{L} h_k \\tau \\right) \\right) \\right \\Vert\n",
    "+\n",
    "\\left \\Vert \\mathcal{R}_1 \\left( \\prod_{k=1}^{L} \\exp (-i h_k \\tau) \\right) \\right \\Vert\n",
    "\\\\\n",
    "\\leq~ &\n",
    "\\left \\Vert \\mathcal{R}_1 \\left( \\exp \\left( \\vert \\tau \\vert \\cdot \\left \\Vert \\sum_{k=1}^{L} h_k \\right \\Vert \\right) \\right) \\right \\Vert\n",
    "+ \n",
    "\\left \\Vert \\mathcal{R}_1 \\left( \\exp \\sum_{k=1}^{L} \\left( \\vert \\tau \\vert \\cdot \\Vert h_k \\Vert \\right) \\right) \\right \\Vert\n",
    "\\\\\n",
    "\\leq~&\n",
    "2 \\mathcal{R}_1 \\left( \\exp ( \\vert \\tau \\vert L \\Lambda ) \\right )\n",
    "\\\\\n",
    "\\leq~&\n",
    " ( \\vert \\tau \\vert L \\Lambda )^2 \\exp ( \\vert \\tau \\vert L \\Lambda ),\n",
    "\\end{aligned}\n",
    "\\tag{11}\n",
    "$$\n",
    "\n",
    "其中 $\\Lambda = \\max_k \\Vert h_k \\Vert$。随后考虑完整的演化时间 $t = r \\cdot \\tau$，那么模拟长度为 $t$ 的时间演化算符时的误差为：\n",
    "\n",
    "$$\n",
    "\\begin{aligned}\n",
    "\\left \\Vert \\left ( \\exp\\left(-i\\sum_{k=1}^L h_k \\tau \\right)\\right)^r - \\left (\\prod_{k=1}^{L} \\exp\\left(-i h_k \\tau \\right) \\right)^r \\right \\Vert\n",
    "\\leq ~&\n",
    "r \\left \\Vert \\exp\\left(-i\\sum_{k=1}^L h_k \\tau\\right) - \\prod_{k=1}^{L} \\exp\\left(-i h_k \\tau \\right) \\right \\Vert\n",
    "\\\\\n",
    "=~& r (  \\tau L \\Lambda )^2 \\exp ( \\vert \\tau \\vert L \\Lambda )\n",
    "\\\\\n",
    "=~& \\frac{(  t L \\Lambda )^2}{r} \\exp \\left( \\frac{\\vert t \\vert L \\Lambda}{r} \\right).\n",
    "\\end{aligned}\n",
    "\\tag{12}\n",
    "$$\n",
    "\n",
    "这里用到了量子电路中误差线性累积的结论，即 $\\Vert U^r - V^r \\Vert \\leq r\\Vert U - V \\Vert$，不熟悉这一结论的读者可以参考 [4] 中的 4.5.3 节。至此，我们就计算出了 product formula 对于一段完整的演化时间 $t$ 的模拟误差上界，即 (7) 式中的二阶项 $O(t^2/r)$。 \n",
    "\n",
    "实际上，我们还可以通过 Suzuki 分解的方法来进一步优化对每一个“时间块”内的时间演化算符 $e^{-iH\\tau}$ 的模拟精度。对于哈密顿量 $H = \\sum_{k=1}^{L} a_k h_k$，其时间演化算符的 Suzuki 分解可以写为\n",
    "\n",
    "$$\n",
    "\\begin{aligned}\n",
    "S_1(\\tau) &= \\prod_{k=0}^L \\exp ( -i h_k \\tau)，\n",
    "\\\\\n",
    "S_2(\\tau) &= \\prod_{k=0}^L \\exp ( -i h_k \\frac{\\tau}{2})\\prod_{k=L}^0 \\exp ( -i h_k \\frac{\\tau}{2})，\n",
    "\\\\\n",
    "S_{2k}(\\tau) &= [S_{2k - 2}(p_k\\tau)]^2S_{2k - 2}\\left( (1-4p_k)\\tau\\right)[S_{2k - 2}(p_k\\tau)]^2,\n",
    "\\end{aligned}\n",
    "\\tag{13}\n",
    "$$\n",
    "\n",
    "对于 $k > 1, k\\in\\mathbb{Z}$，其中 $p_k = 1/(4 - 4^{1/(2k - 1)})$。先前推导的 product formula 实际上只使用了一阶的 Suzuki 分解 $S_1(\\tau)$ 来对每个“时间块”进行模拟，因此也被称为一阶 Suzuki product formula，或者简称一阶 product formula。Suzuki product formula 在物理中也经常被称为 Trotter-Suzuki 分解方法。对于更高阶的 product formula 而言，使用 (10-12) 中类似的计算，可以证明第 $2k$ 阶 product formula 的整体误差上界为：\n",
    "\n",
    "$$\n",
    "\\epsilon\\left(e^{-iHt}, \\left(S_{2k}(\\tau)\\right)^r\\right)\n",
    "\\leq\n",
    "\\frac{(2L5^{k-1}\\Lambda\\vert t \\vert)^{2k+1}}{3r^{2k}} \\exp \\left( \\frac{2L5^{k-1}\\Lambda\\vert t \\vert}{r} \\right),\n",
    "~ k > 1.\n",
    "\\tag{14}\n",
    "$$\n",
    "\n",
    "在得到了模拟误差上界的基础上，还可以进一步计算达到一定最小精度 $\\epsilon$ 时所需要的电路深度的下界。需要指出的是，(12) 和 (14) 中给出的误差上界的计算是十分宽松的。近年来，许多工作都进一步地给出了更加紧致的上界 [3, 5-6]。此外，也有人提出了不基于 Suzuki 分解的 product formula [7]。"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f3bf6fb2",
   "metadata": {},
   "source": [
    "![image.png](./figures/trotter_suzuki_circuit.png)\n",
    "<div style=\"text-align:center\">图1：Suzuki product formula 电路的示意图</div>"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e0c26758",
   "metadata": {},
   "source": [
    "### 利用 Paddle Quantum 验证基于 Suzuki-product formula 的时间演化电路\n",
    "\n",
    "尽管人们不断的在优化 Suzuki-product formula 的误差上界，但是在实际使用中，真实的误差往往和理论上给出的上界有一定差距，也就是说，我们现在能计算的理论上的 product formula 误差依然只是一个十分宽松的上界 [3]。因此，对于一个真实的物理系统，我们往往需要通过数值实验来计算其真实的误差，从而给出一个经验性的误差上界（empirical bound）。这样的数值试验对于计算在一定精度下模拟某一时间演化过程所需要的电路深度十分重要。\n",
    "\n",
    "在刚刚介绍的 `construct_trotter_circuit` 函数中，用户也可以通过改变传入参数 `tau`、`steps`、`order` 来实现对任意阶 Suzuki product formula 电路的创建。下面，我们就利用该函数的这一功能来验证一下 Suzuki product formula 的误差。\n",
    "\n",
    "继续使用先前的哈密顿量："
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "86de7992",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "H = 1.0 Z0, Z1\n",
      "1.0 X0\n",
      "1.0 X1\n"
     ]
    }
   ],
   "source": [
    "print('H =', h_2)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2e183412",
   "metadata": {},
   "source": [
    "这里我们通过改变 `tau` 和 `steps` 两个参数来将 $t=1$ 的演化过程进行拆分。（提示：$\\tau \\cdot n_{\\rm steps} = t$）"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "240ed09f",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "--*-----------------*----Rx(0.667)----*-----------------*----Rx(0.667)----*-----------------*----Rx(0.667)--\n",
      "  |                 |                 |                 |                 |                 |               \n",
      "--x----Rz(0.667)----x----Rx(0.667)----x----Rz(0.667)----x----Rx(0.667)----x----Rz(0.667)----x----Rx(0.667)--\n",
      "                                                                                                            \n",
      "电路的酉算符和时间演化算符之间的保真度为： 0.984\n"
     ]
    }
   ],
   "source": [
    "# 将长度为 t 的时间演化过程拆为 r 份\n",
    "r = 3\n",
    "t = 1\n",
    "cir = Circuit()\n",
    "# 搭建时间演化电路，tau 为每个“时间块”的演化时间，故为 t/r，steps 是“时间块”的重复次数 r\n",
    "construct_trotter_circuit(cir, h_2, tau=t/r, steps=r)\n",
    "print(cir)\n",
    "print('电路的酉算符和时间演化算符之间的保真度为： %.3f' \n",
    "      % gate_fidelity(cir.unitary_matrix().numpy(), linalg.expm(-1 * 1j * h_2.construct_h_matrix())))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4708cad9",
   "metadata": {},
   "source": [
    "我们发现，通过将 $t=1$ 的模拟时间拆分为三个“时间块”，模拟误差被成功减小了。\n",
    "\n",
    "如果进一步地增加拆分后的“时间块”的数量的话，误差还可以被进一步减小："
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "abc3742d",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZ4AAAEJCAYAAACkH0H0AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAkmElEQVR4nO3deZRcZZ3/8fe3qveks5B9gwQSAiEqEUTFAfkZHBaDMC5AFEaEEcVxZpiDOIAL4wCuOIMelREHBpUxhF1h2BRGOCJigERJZSEhLOmu7EtVd6c73V31/f1Rt5Pu6upOd9Jdt5bP65w6XXXvrXu/VenUp5/nPvVcc3dERETyJRJ2ASIiUl4UPCIiklcKHhERySsFj4iI5JWCR0RE8qoi7AKKwfjx433mzJlhlyEiUlReeuml7e4+IXu5gmcAZs6cyYsvvhh2GSIiRcXM3sy1XF1tIiKSVwoeERHJKwWPiIjklYJHRETySsEjIiJ5VXaj2sxsBPBjoB34nbv/T8gliYiUlZIIHjO7A1gEbHX3+d2Wnwl8H4gC/+Xu3wI+Atzn7g+b2VJAwSMiPLS8ke8+sZb47lamjqnl6jPmct6CaWVXQz7qsFK4LIKZnQo0Az/vCh4ziwKvAh8EGoBlwGLgXOAxd19hZr90908caP/19fV+wgkn9Fh2/vnn8/nPf549e/Zw9tln93rOJZdcwiWXXML27dv52Mc+1mv9FVdcwQUXXMDGjRu5+OKLe62/6qqrOOecc1i7di2f/exne63/yle+wumnn86KFSu48sore63/xje+wcknn8wf/vAHrrvuul7rb7nlFo4//nh++9vfcuONN/Za/5Of/IS5c+fy8MMP873vfa/X+l/84hfMmDGDpUuXcuutt/Zaf9999zF+/HjuvPNO7rzzzl7rH330Uerq6vjxj3/MPffc02v97373OwBuvvlmHnnkkR7ramtreeyxxwC44YYbeOqpp3qsHzduHPfffz8A1157Lc8//3yP9dOnT+euu+4C4Morr2TFihU91h999NHcdtttAFx++eW8+uqrPdYff/zx3HLLLQBcdNFFNDQ09Fj/3ve+l29+85sAfPSjH2XHjh091i9cuJCvfvWrAJx11lm0trb2WL9o0SK++MUvAnDaaaeRbah+925/6i9845GVpKpGEW1PMvatZxm5Y03ef/eaxx3DrsNP3VfHFz84hyvOfldef/eaxx3DjiPPxKOV+7arIM2Ydf/LyB1r9i0bzt+95nHHsHP2WaRtf3vAUh2M2/D4vhry8bs3e+Firn3gFVo7Ur3quHTh2wf1u/fMM8+85O4nZm9XEi0ed3/WzGZmLT4JWO/uGwDM7G4yodMATAdW0M85LjO7HLgcoLq6euiLlrK1vmMs7/vW08R3txI99m8Z8+YzPT7c8uGh5Y185/8aSFWPBiBVPZodR545ZPt3d9LuuEVxi4AZTgQsws7WFPHdrWxrdXZPPoHE4afgkcp9dfz7c9vYW/UqkSS0jZyWea4ZYGDGC282sb5lK6sTUfaMnY0HyyGz3WOrtjNiZBvLd1fTNGE+WGTfOjB+8UID0YoK/pSsJzHlJNwgMfXdPUIHoJMIO2b9Ne11E8EAjIqKSm54ZBVpd17YO50dMxcSrMTNaK6p4Zr7/0LanWU+h21Hjdm3HjNa6ur4wi9fxh2WV76d3XOOCN4baBs9E7eeH8kerWT7kWfRPPHtAPy+sp4Lb3sed1g9YSF7R7Xv2zfAk17P6h89hwPrZpxD55RUj/UPttbz3C3P4g5vHL2YdFfDI6jhFztH0nzvn0mlezZIPFrJrsNPBXYfzK9DLyXR4gEIgueRbi2ejwFnuvvfBY8vBt4N/AvwQ6AN+P1AzvGceOKJrpkLZCg8tLyx11+TtZVRvvmRt+3rykinnY50mo6U096ZpiOVpr0zTXtq//3MT88s63q8b53T3pnK/Oyxfdd2zkPLG3vU0KW6IsLxM8aQSjud6Ux4dKY8eJwm7dCZTpNKdVufdlIpJ9V1P7iViqpoBAuyLWJGxAwjeBzJ3I+YBdsYEQMj+Bksj2T97HpO9+et3pTss4YTjhi775hGZgc9j5tZHuRLj+NYt3p73Kf3c82Mh/8cz1mDAa9/60ODeu/MrHRbPH2wHMvc3VuAT+e7GAnfUPVbp9POno4ULXs7adnbyZ721P6f7Z3B8hR72jtpaU+xZ2/mZ8veTp5es5W9neke+2vtSPHP96zgugdfoSOVCY6hVhk1KqMRKqMRqioiOUMH2FdbdWWEWjMqIkY0Egl+Zm4VESMSsX3Lej6O9Fof7WfbL9775z5r/tmlJxHp9mEfMYhGbN+HezRi+9dFIGo51kWCx1nrzCzYJrP/D9z8O+KJtl41TBtTy3PXfGBo/hEO4H3feprG3a29lk8bU8v9V5yclxoAXn5zV846po6pHbJjlHLwNAAzuj2eDuSOcil5Dy1v5JoH/kJbR+aDtXF3K1ff92ee37CduZNGsae9k+ausMgKjeYgVLrW9fWhnUtl1BhRXcGIqgrqqqK9QqeLO3zipMOpqtgfDlXRCJVRo6oiGvzsWhahMrhfVWH7tq+Mdi2LdFtmVEYiRCI9/w7r70Nu6WffO4h39tD8x29e7bOO9x/da27JYfOlM4/J2RK9+oy5eavh6jPmhl5Dvuoo5eBZBswxs1lAI3AhcMCBBFKc2jpSbE60EU+0sjnRxqZEG5sSrWzanbm/ZnOS7N6fjpSzdNn+k7NV0Qh11dF9ITGiuoIR1VHGjqhjRFWUuuoKRnQtr6rIsW1Fj+3qqiqoquh5GrG/D/yvLJo3LO9NLuX0ITcQXS3fMEeUFUIN+aqjJM7xmNkS4DRgPLAFuN7dbzezs4FbyAynvsPdbzqY/escz6E7lG6uvkJlc6KN+O42Nifb2NnS3ut5Y+sqmTy6limja3h6zdac+zZgxdf+mtqqaK+QGA4DOceTL+UydFfC09c5npIInuGm4Dk0/X3Ynjl/cs8WygBDZUxdJVOCUNl/Cx6PqWXyqBpqq6L7tu+vpZGvPvwu+qCVcqHgOQQKnkPT14d+xOjV/QWZUJk8qoapY2qZPLqGqaNrmDy6NviZCZjuoTIQhdTSECkX5TiqTULW3pnmqdVbcoYOZELnqg8ezZQx+1suk0fXUFc19L+WhdJ/LiIKHhkG67c2c8+LG7n/pQZ2tLT32bKZNqaWf1g4J291nbdgmoJGpAAoeGRItLan+N9XNrF02Vsse2MXFRHj9GMnccFJM9jV3M6XH1oZ+sglESkMCh45JCsbE9y97C1+tTxO095OZo0fwTVnHcNH3jmNifU1+7aLREzdXCICKHjkICRaO/j1ikbuXraRWDxJdUWED71tChe8awYnzToMs96TRqibS0S6KHhkQNydP72+k6XLNvK/r2xib2eaeVNGccO5x/Hh46cxurbywDsREUHBIwewrWkvD7zcwNJlG9mwvYX66go+dsJ0Fp90OPOnjQ67PBEpQgoe6SWVdp5dt42lf9rIb1dvoTPtvGvmWD7//2Zz9tsmD8twZxEpH/oEkX0adu3hnhcbuPfFjWxKtDFuRBWX/tUszj9xBrMnjgy7PBEpEQqeMtfemeY3q7Zw97K3+P367QCcMmcCX100j9OPnZSX+ctEpLwoeMrU+q1NLF22kftfbmRnSztTR9fwjx+Yw8dPnM70sXVhlyciJUzBU+K6T0g5eXQNp82dwLotzbz4ZuZLnh+cN4kL3jWDU+ZMIBrJde08EZGhpeApYdkTY25KtLHkTxuZUF/FdWcfw98smM6E+uqQqxSRcqPgKWHffWJtzqtlVkUjXH7qUSFUJCICOnNcwuJ9zAod39372vIiIvmi4ClhU8fUDmq5iEg+KHhK2NVnzO01YECzQotI2BQ8Jey8BdOYMbaWqmgEI3P9G11xU0TCpsEFJSyVdrY27WXxSTP4+rnzwy5HRARQi6ekvbGjhT3tKY6bqsk8RaRwKHhKWCyeBOC4aaNCrkREZD8FTwmLxRNURo05E+vDLkVEZB8FTwlbFU9y9KR6TfQpIgVFn0glyt2JxZMcN1XdbCJSWBQ8JWpzso2dLe0aWCAiBUfBU6JWNgYDC9TiEZECo+ApUbF4AjM4doqCR0QKi4KnRMXiSWaNG8GIan1HWEQKi4KnRK2KJ5mnbjYRKUAKnhK0q6Wdxt2tGlggIgVJwVOCVm3SwAIRKVwKnhIUiycABY+IFCYFTwmKxZNMHlXDuJHVYZciItKLgqcExeJJ5mtiUBEpUAqeEtPanmLDtmbmaWCBiBQoBU+JWb05Sdp1fkdECpeCp8TsuwaPgkdECpSCp8SsiicYXVvJtDG1YZciIpKTgqfEdF0KwczCLkVEJCcFTwnpSKVZs6lJ3WwiUtAUPCVk/dZm2lNpTZUjIgVNwVNCNLBARIqBgqeExOIJaiojHDlhZNiliIj0qWyDx8yONLPbzey+sGsZKrF4kmMmjyIa0cACESlceQseM/snM1tpZjEzu/IQ9nOHmW01s5U51p1pZmvNbL2ZXdPfftx9g7tfdrB1FJp02lkdjGgTESlkeQkeM5sPfAY4CXgHsMjM5mRtM9HM6rOWzc6xuzuBM3McIwr8CDgLmAcsNrN5ZvY2M3sk6zZxSF5YAdm4aw9Nezs1sEBECl6+WjzHAn909z3u3gk8A/xN1jbvB35lZjUAZvYZ4AfZO3L3Z4GdOY5xErA+aMm0A3cD57r7K+6+KOu2dSBFm9k5ZnZbIpEY8AsNS9fAAk0OKiKFLl/BsxI41czGmVkdcDYwo/sG7n4v8Dhwt5l9ErgUOH8Qx5gGbOz2uCFYllNQy38CC8zs2lzbuPvD7n756NGF34qIxRNEI8bRk+oPvLGISIgq8nEQd19tZt8GfgM0A38GOnNs9x0zuxu4FTjK3ZsHcZhcZ9S9n5p2AJ8bxP4LWiyeZM7EkdRURsMuRUSkX3kbXODut7v7O939VDJdZeuytzGzU4D5wIPA9YM8RAM9W1HTgfhBllt0YvEk8zSwQESKQD5HtU0Mfh4OfARYkrV+AfBT4Fzg08BhZnbjIA6xDJhjZrPMrAq4EPj1UNRe6LY2tbGtaa8GFohIUcjn93juN7NVwMPA37v7rqz1dcDH3f01d08DnwLezN6JmS0BngfmmlmDmV0GEAxa+ALwBLAauMfdY8P3cgqHZiwQkWKSl3M8AO5+ygHWP5f1uINMCyh7u8X97ONR4NGDrbFYrQqCR11tIlIMynbmglKysjHB4YfVMaqmMuxSREQOSMFTAmKasUBEioiCp8gl2zp4a+ceBY+IFA0FT5FbtW9ggUa0iUhxUPAUOY1oE5Fio+ApcrF4gvEjq5k4qibsUkREBkTBU+RWxZOaGFREioqCp4i1daRYt7VZ3WwiUlQUPEXs1S1NpNKugQUiUlQUPEVMAwtEpBgpeIpYLJ6gvrqCGWPrwi5FRGTAFDxFLBZPcuzUUUQiuS5FJCJSmBQ8RSqVdtZsalI3m4gUHQVPkXp9ezOtHSkNLBCRoqPgKVIrGzWwQESKk4KnSMXiCaoqIsyeODLsUkREBkXBU6Ri8SRzJ9VTGdU/oYgUF31qFSF31zV4RKRoKXiKUOPuVhKtHQoeESlKCp4itG/Ggmka0SYixUfBU4Ri8SQRg2Mnq8UjIsVnQMFjZg+Y2XlmVjncBcmBrYonOHLCSGqromGXIiIyaANt8TwHfA3YbGa3mtnJw1iTHIAGFohIMRtQ8Lj799z9ncCpwG5giZmtN7OvmdlRw1mg9LSzpZ1NiTYFj4gUrUGd43H3mLtfC1wEtADXAy+b2W/N7B3DUaD0FIsnADRVjogUrQEHj5nNNbMbzOw14DZgKTATmAQ8Cjw0HAVKT7oGj4gUu4qBbGRmL5IJmaXAJ9z9haxN/t3M/mGIa5McYvEk08bUMqauKuxSREQOyoCCB/gW8Gt3b+9rA3efNTQlSX9i8QTz1NoRkSI20K62L+cKnaAlJHnSsreT17e3qJtNRIraQIOn18g1MzPgyKEtR/qzelMSdw0sEJHi1m9Xm5n9PLhb3e1+l5lAbDiKktw0sEBESsGBzvG81sd9J/Ol0nuHvCLpUyyeYGxdJVNG14RdiojIQes3eNz96wBm9kd3fyI/JUlfYvEk86eNJtPLKSJSnPoMHjM71d2fDR52mNkHcm3n7k8PS2XSQ3tnmle3NHHpX2nwoIgUt/5aPD8G5gf3b+9jG0cDDPJi3dYmOlKugQUiUvT6DB53n9/tvv7MDpkGFohIqdD1eIrEqniSuqoos8aNCLsUEZFD0t85no1kutL65e6HD2lFklMsnuDYKaOIRDSwQESKW3/neC7KWxXSr3TaWRVP8tETpoddiojIIevvHM8z+SxE+vbmzj20tKd0fkdESsJAL31dbWY3mdkGM0sEy/7azL4wvOUJ6Bo8IlJaBjq44D/IDK3+JPvP+8SAK4ajKOkpFk9SETHmTBoZdikiIodsoJdF+Btgtru3mFkawN0bzWza8JUmXVY2JpgzqZ7qimjYpYiIHLKBtnjayQopM5sA7BjyiqQH98zAAp3fEZFSMdDguRf4mZnNAjCzKcAPgbuHqzDJ2JLcy46WdgWPiJSMgQbPdcAbwCvAGGAdEAe+PixVyT5dAwvmT9PAAhEpDQM6xxNcffRK4Mqgi227ux/wy6Vy6GLxJGZw7BS1eESkNPQ3c0F/k3/Wd03N7+4bhroo2S8WTzBz3AhGVg90HIiISGHr79NsPZmh08b+IdRd87V0b+1oqNUwisWTvGPGmLDLEBEZMn2e43H3iLtH3T0C/B2ZgQRzgRrgGOCXwGV5qbJMJfZ00LCrVQMLRKSkDLT/5gZgjru3Bo/XmdlngVeBO4ejMIHYJs1YICKlZ6Cj2iLAzKxlR1DE3WxmdqSZ3W5m94VdS19W6Ro8IlKCBjNlztNm9g0zu8LMvgE8FSwfEDP7ZzOLmdlKM1tiZjUHU7CZ3WFmW81sZY51Z5rZWjNbb2bX9Lcfd9/g7gXdVRiLJ5k0qprxI6vDLkVEZMgMKHjc/bvAp4FJwIeBycCl7v6dgTw/mFrnH4ETgyubRoELs7aZaGb1Wctm59jdncCZOY4RBX4EnAXMAxab2Twze5uZPZJ1mziQusMWiyfUzSYiJWfAY3Td/XHg8UM8Vq2ZdQB1ZL6A2t37gSvM7Gx3bzOzz5CZI+7srDqeNbOZOfZ/ErC+a3i3md0NnOvu3wQWHUzBZnYOcM7s2bnyb3i1daR4bVsLZxw3Oe/HFhEZTv19j+fL7n5TcP/f+trO3b92oIMEE4reDLwFtAJPuvuTWdvcG0zJc7eZ3QtcCnxwYC8DgGnAxm6PG4B397WxmY0DbgIWmNm1QUBl1/0w8PCJJ574mUHUMSTWbG4ilXad3xGRktNfi+frZD6YAY4iM1HoQTGzscC5wCxgN3CvmV3k7nd1387dvxO0VG4FjnL35sEcJseyPmdXcPcdwOcGsf+8WtmoEW0iUpr6C5493e6f4+6H8qf36cDr7r4NwMweAE4GegSPmZ1C5ro/DwLXA4O50FwDMKPb4+n07s4rGrF4klE1FUwfWxt2KSIiQ6rfmQvM7HtkLvhWYWafJkerwt3vGMBx3gLeY2Z1ZLraFgIvdt/AzBYAPwU+BLwO3GVmN7r7Vwb0SmAZMCformskM3jhEwN8bsFZFQws6JqaSESkVPQXPBcCXwIWA1XA3+bYxoEDBo+7vxB8X+ZloBNYDtyWtVkd8HF3fw3AzD4FXJK9LzNbApwGjDezBuB6d7/d3TuDS3E/QWbU3B3uHjtQbYWoM5VmzeYmLn7PEWGXIiIy5PoMHnd/lcxUOZjZU+6+8FAO5O7Xk+k+62v9c1mPO8i0gLK3W9zPPh4FHj2EMgvCa9ta2NuZ5rhpGlggIqVnoN/jOaTQkcHpugaPBhaISCka6MwFkkexeJLqighHjh8RdikiIkNOwVOAYvEEx0wZRUVU/zwiUnr0yVZg3J1V8aS+OCoiJUvBU2AadrWSbOtU8IhIyVLwFBgNLBCRUqfgKTCxeJJoxDhmcv2BNxYRKUIKngITiyc5asIIaiqL9hp7IiL9UvAUmJWNugaPiJQ2BU8B2da0l61NezWwQERKmoKngGhggYiUAwVPAYnFkwDMU4tHREqYgqeArIonmXFYLaNrK8MuRURk2Ch4CkgsnuC4KepmE5HSpuApEE1tHbyxY48GFohIyVPwFIjVm5oAdA0eESl5Cp4CoRFtIlIuFDwFIhZPMn5kFRPrq8MuRURkWCl4CkQsnmTe1NGYWdiliIgMKwVPAdjbmWLdliYNLBCRsqDgKQDrtjTTmXYFj4iUBQVPAdDAAhEpJwqeArCyMcnI6gqOOKwu7FJERIadgqcAxOIJ5k0ZRSSigQUiUvoUPCFLpZ3Vm5o0MaiIlA0FT8he395Ca0dKAwtEpGwoeEKmgQUiUm4UPCFbFU9SFY0wZ9LIsEsREckLBU/IYvEkR08eSWVU/xQiUh70aRcid9c1eESk7Ch4QrQp0cauPR26FIKIlBUFT4hi8SSARrSJSFlR8IQoFk9gBsdMVvCISPlQ8IQoFk8ya/wIRlRXhF2KiEjeKHhCtCqe1Pd3RKTsKHhCsqulncbdrTq/IyJlR8ETkq6BBfPV4hGRMqPgCcn+qXLU4hGR8qLgCUksnmTq6BrGjqgKuxQRkbxS8IQkFk8wT91sIlKGFDwh2NPeyYbtLepmE5GypOAJwepNTbjr/I6IlCcFTwhWdQ0smKauNhEpPwqeEMTiScbUVTJ1dE3YpYiI5J2CJwSxeJLjpo7CzMIuRUQk7xQ8edaRSrN2c5OmyhGRsqXgybP1W5tpT6U1sEBEypaCJ890DR4RKXcKnjyLxRPUVkaZNX5k2KWIiIRCwZNnsXiSY6fUE41oYIGIlCcFTx6l065r8IhI2VPw5NFbO/fQvLdT53dEpKwpePJo/8ACtXhEpHwpePIoFk9QETGOnqyBBSJSvhQ8eRSLJ5k9cSTVFdGwSxERCY2CJ49iGlggIqLgyZetyTa2N+/VwAIRKXsKnjzRjAUiIhllGzxmdqSZ3W5m9+XjeLHgGjzzFDwiUubyEjxmNtfMVnS7Jc3syoPc1x1mttXMVuZYd6aZrTWz9WZ2TX/7cfcN7n7ZwdRwMGLxJEeMq6O+pjJfhxQRKUgV+TiIu68FjgcwsyjQCDzYfRszmwi0untTt2Wz3X191u7uBH4I/Dzr+VHgR8AHgQZgmZn9GogC38zax6XuvvXQXtXgxOJJ5k9Ta0dEJIyutoXAa+7+Ztby9wO/MrMaADP7DPCD7Ce7+7PAzhz7PQlYH7Rk2oG7gXPd/RV3X5R1y2voJNs6eGvnHo1oExEhnOC5EFiSvdDd7wUeB+42s08ClwLnD2K/04CN3R43BMtyMrNxZvafwAIzu7aPbc4xs9sSicQgyuhtlQYWiIjsk9fgMbMq4MPAvbnWu/t3gDbgVuDD7t48mN3n2mVfG7v7Dnf/nLsf5e7ZXXFd2zzs7pePHn1oLRVNlSMisl++WzxnAS+7+5ZcK83sFGA+mfM/1w9y3w3AjG6PpwPxgylyqMUaE0ysr2ZCfXXYpYiIhC7fwbOYHN1sAGa2APgpcC7waeAwM7txEPteBswxs1lBy+pC4NeHWO+QyMxYoG42ERHIY/CYWR2ZEWcP9LFJHfBxd3/N3dPAp4DsAQiY2RLgeWCumTWY2WUA7t4JfAF4AlgN3OPusaF/JYPT1pFi/bZmdbOJiATyMpwawN33AOP6Wf9c1uMOMi2g7O0W97OPR4FHD6HMIbd2cxOptKvFIyISKNuZC/JFAwtERHpS8AyzWDxBfU0FMw6rDbsUEZGCoOAZZrF4knlTRmGWa7S3iEj5UfAMo1TaWbNZ1+AREelOwTOMNmxrpq0jrYEFIiLdKHiGyUPLGzn/J88D8O3H1/DQ8saQKxIRKQx5G05dTh5a3si1D7xCa0cKgK1Ne7n2gVcAOG9Bn9PHiYiUBbV4hsF3n1i7L3S6tHak+O4Ta0OqSESkcCh4hkF8d+uglouIlBMFzzCYOib3d3b6Wi4iUk4UPMPg6jPmUlsZ7bGstjLK1WfMDakiEZHCocEFw6BrAMF3n1hLfHcrU8fUcvUZczWwQEQEBc+wOW/BNAWNiEgO6moTEZG8UvCIiEheKXhERCSvFDwiIpJXCh4REckrc/ewayh4ZrYNeDPsOg7ReGB72EUUCL0XPen96Envx36H+l4c4e4TshcqeMqEmb3o7ieGXUch0HvRk96PnvR+7Ddc74W62kREJK8UPCIiklcKnvJxW9gFFBC9Fz3p/ehJ78d+w/Je6ByPiIjklVo8IiKSVwoeERHJKwVPCTOzGWb2f2a22sxiZvZPYddUCMwsambLzeyRsGsJm5mNMbP7zGxN8Hvy3rBrCouZ/XPw/2SlmS0xs5qwa8onM7vDzLaa2cpuyw4zs9+Y2brg59ihOJaCp7R1Ale5+7HAe4C/N7N5IddUCP4JWB12EQXi+8Dj7n4M8A7K9H0xs2nAPwInuvt8IApcGG5VeXcncGbWsmuAp9x9DvBU8PiQKXhKmLtvcveXg/tNZD5UyvoiQWY2HfgQ8F9h1xI2MxsFnArcDuDu7e6+O9SiwlUB1JpZBVAHxEOuJ6/c/VlgZ9bic4GfBfd/Bpw3FMdS8JQJM5sJLABeCLmUsN0CfAlIh1xHITgS2Ab8d9D1+F9mNiLsosLg7o3AzcBbwCYg4e5PhltVQZjk7psg84csMHEodqrgKQNmNhK4H7jS3ZNh1xMWM1sEbHX3l8KupUBUAO8EbnX3BUALQ9SVUmyCcxfnArOAqcAIM7so3KpKl4KnxJlZJZnQ+R93fyDsekL2PuDDZvYGcDfwATO7K9ySQtUANLh7Vyv4PjJBVI5OB153923u3gE8AJwcck2FYIuZTQEIfm4dip0qeEqYmRmZ/vvV7v7vYdcTNne/1t2nu/tMMieOn3b3sv2r1t03AxvNbG6waCGwKsSSwvQW8B4zqwv+3yykTAdaZPk18Kng/qeAXw3FTiuGYidSsN4HXAy8YmYrgmXXufuj4ZUkBeYfgP8xsypgA/DpkOsJhbu/YGb3AS+TGQ26nDKbOsfMlgCnAePNrAG4HvgWcI+ZXUYmnD8+JMfSlDkiIpJP6moTEZG8UvCIiEheKXhERCSvFDwiIpJXCh4REckrBY9InpjZG2Z2ekjHnmRmz5pZk5l9L4waRLroezwi5eFyYDswygfxHQozOw24y92nD1NdUobU4hEpMsHsyYN1BLBqMKEjMlwUPFLWgu6vL5rZX8wsYWZLuy4AZmaXmNnvs7Z3M5sd3L/TzH5sZo+ZWbOZPWdmk83sFjPbFVxcbUHWId9lZquC9f/d/WJjZrbIzFaY2W4z+4OZvT2rzn8xs78ALbnCx8xONrNlwetYZmYnd9VJZrqTLwV19uruM7Ozg7qazKwxeE9GAI8BU4PnNZvZVDOLmNk1Zvaame0ws3vM7LBgPzOD9+hyM4ub2SYzu6rbcU4ysxfNLGlmW8ys7KdyKkvurptuZXsD3gD+RGZG4sPIzM/1uWDdJcDvs7Z3YHZw/04y3VcnADXA08DrwN+SuZDYjcD/ZR1rJTAjONZzwI3BuneSmYDx3cFzPxVsX93tuSuC59bmeB2HAbvITJFUASwOHo/rVuuN/bwPm4BTgvtjgXcG908jM5Fo922vBP4ITAeqgZ8AS4J1M4P3aAkwAngbmUsvnB6sfx64OLg/EnhP2L8DuuX/phaPCPzA3ePuvhN4GDh+EM990N1fcvc24EGgzd1/7u4pYCmZayB190N33xgc6yYyAQHwGeAn7v6Cu6fc/WfAXjJXju1e50Z3b81Rx4eAde7+C3fvdPclwBrgnAG+jg5gnpmNcvddHlxAsA+fBb7s7g3uvhf4V+BjWa2wr7t7i7u/Avx3t9fZAcw2s/Hu3uzufxxgfVJCFDwisLnb/T1k/hIfqC3d7rfmeJy9r43d7r9JpqUFmXMwVwXdbLvNbDeZ1s3UPp6bbWqwv+7eZOBXnP0ocDbwppk9Y2bv7WfbI4AHu9W5GkgBk/qotfvrvAw4GlgTdAcuGmB9UkIUPCJ9ayFzCWQAzGzyEOxzRrf7h7P/8sobgZvcfUy3W13QcunS38CAOJlA6O5woHEgRbn7Mnc/l8wVJh8C7unnmBuBs7JqrfHMVTy75Hyd7r7O3RcHx/k2cF+5XvW0nCl4RPr2Z+A4Mzs+GATwr0Owz783s+nByfjryHTHAfwU+JyZvdsyRpjZh8ysfoD7fRQ42sw+YWYVZnYBMA945EBPNLMqM/ukmY32zEXQkmRaMJBpwY0zs9HdnvKfwE1mdkTw/Almdm7Wbr8aXNvmODKXWlgabHuRmU1w9zSwO9g2hZQVBY9IH9z9VeDfgN8C64Df9/+MAfkl8CSZa99sIDMAAXd/kcx5nh+SGRSwnszghoHWugNYBFwF7AC+BCxy9+0D3MXFwBtmlgQ+B1wU7HcNmYECG4KutanA98lcIOxJM2siM9Dg3Vn7eyZ4DU8BN7v7k8HyM4GYmTUH+7kwOD8mZUTX4xGRIWNmM8mM7Kt0986Qy5ECpRaPiIjklYJHRETySl1tIiKSV2rxiIhIXil4REQkrxQ8IiKSVwoeERHJKwWPiIjk1f8Ha8R4cDtyy6QAAAAASUVORK5CYII=",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# 导入一些作图和计算误差需要的额外包\n",
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "\n",
    "def get_fid(n_steps):\n",
    "    t = 1\n",
    "    cir = Circuit()\n",
    "    construct_trotter_circuit(cir, h_2, tau=t/n_steps, steps=n_steps)\n",
    "    return gate_fidelity(cir.unitary_matrix(), linalg.expm(-1 * 1j * h_2.construct_h_matrix()))\n",
    "plt.axhline(1, ls='--', color='black')\n",
    "plt.semilogy(np.arange(1, 11), [get_fid(r) for r in np.arange(1, 11)], 'o-')\n",
    "plt.xlabel('number of steps', fontsize=12)\n",
    "plt.ylabel('fidelity', fontsize=12)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f8381126",
   "metadata": {},
   "source": [
    "不仅如此，我们还可以通过改变 product formula 的阶数来减小模拟误差。目前，`construct_trotter_circuit` 函数支持通过参数 `order` 来实现任意阶数的 Suzuki product formula。下面，就让我们来分别计算一下一阶和二阶时间演化电路的误差随着 $\\tau$ 大小的变化，并和上文中计算的理论误差上界进行对比："
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "id": "d3256b18",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZEAAAENCAYAAADOhVhvAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAzFUlEQVR4nO3dd3iVVbr38e9KSEKvoSWhCoROggEcUYrSLDQ7drGOZ5xxZvSoc+aMTvHoHD2+jmMZsY6OgoqQgA0VC1iGmtCESA/ZgQCBhJa+1/vHk2hAIMnO3nn2zv59riuX7PbkNj7mZq11r3sZay0iIiK+iHA7ABERCV1KIiIi4jMlERER8ZmSiIiI+ExJREREfKYkIiIiPmvkdgD1LTY21nbv3t3tMEREQsaqVav2W2vbn+y1sEkixpjJwORevXqxcuVKt8MREQkZxpidp3otbKazrLULrbW3tWrVyu1QREQajLBJIiIi4n9hk0SMMZONMbMKCgrcDkVEpMEImzURa+1CYGFKSsqtJ75WWlpKdnY2RUVFLkQWPBo3bkxCQgJRUVFuhyIiISJskkjVhfUTZWdn06JFC7p3744xpv6DCwLWWvLy8sjOzqZHjx5uhyMiISJsprNOt7BeVFREu3btwjaBABhjaNeuXdiPxkSkdsImiVQnnBNIJf0MRBqowoPw3YKAXDpsprNERMJSaRHMvho8qyB+NbRK8OvlwyaJnG5NpLZS0z08tiiTnPxC4lo34d6JiUxLjq97kKdRXl5OZGTkKR+fjLUWay0RERpwioQlbznMuwWyvoFLX/J7AoEwms7y12bD1HQPD8xbhye/EAt48gt5YN46UtM9dbruv/71L4YPH05SUhK333475eXlNG/enD/84Q+MGDGCb7/99iePn3jiCQYOHMjAgQN58sknAdixYwf9+vXjzjvvZOjQoezatatOcYlIiLIWPrwPNi6EiY/AoMsC8m3CZiRSU39cuIHvcg6d8vX0rHxKyr3HPVdYWs5/zl3L7OVZJ/1M/7iWPDh5wCmvuXHjRt566y2+/vproqKiuPPOO3njjTc4evQoAwcO5E9/+hPAcY9XrVrFK6+8wrJly7DWMmLECEaPHk2bNm3IzMzklVde4dlnn/XhJyAiDcLS/4MVL8DZd8HP7gzYtwmbJOKv6awTE0h1z9fE4sWLWbVqFcOGDQOgsLCQDh06EBkZyaWXXvrD+6o+/uqrr5g+fTrNmjUD4JJLLmHp0qVMmTKFbt26cdZZZ/kcj4iEuPQ34LM/w6ArYNyfAvqtwiaJnG6zYVWnGzEAjHz0Mzz5hT95Pr51E966/We+xsYNN9zAI488ctzzjz/++HHrHo0bN/7hsbX2lNerTCwiEoY2fwIL7oKeY2DqMxDgNdGwWRPxl3snJtIk6vgF7SZRkdw7MdHna55//vnMnTuXvXv3AnDgwAF27jxl00wARo0aRWpqKseOHePo0aPMnz+fc8891+cYRKQB8KyCt6+Hjv3hitehUXTAv2XYjET8pbIKy5/VWf379+cvf/kLEyZMwOv1EhUVxTPPPHPazwwdOpQbb7yR4cOHA3DLLbeQnJzMjh07fI5DREJY3lZ44wpo1h6ueRcat6yXb2tONy3SEKWkpNgTzxPZuHEj/fr1cymi4KKfhUgIOrIXXhwHJUdg5scQW/etDFUZY1ZZa1NO9lrYTGepi6+INEjFh+GNy+DoPrj6bb8nkOqETRLRoVQi0uCUlThrIHvWw+WvQsJJBwsBpTUREZFQZK1ThbX1M5jyNPSZ6EoYYTMSERFpUD59CNbOgbG/h6HXuRaGkoiISKj59z/g6ychZSaMusfVUJRERERCyfp58NH90PdiuPBxcPkIh7BJIsFenbVjxw4GDhxYb9/voYce4vHHH6+37ycifrB9Kcy/HbqMgEtfhIjTd/KuD2GTRFSdJSIhLXcDzLkG2vSAGbMhqonbEQFhlERCQVlZGTfccAODBw/msssu49ixYyxevJjk5GQGDRrEzJkzKS4uBqB79+7s378fgJUrVzJmzBjAGWHMnDmTMWPG0LNnT5566qkfrv/www+TmJjIuHHjyMzMrPd/PxHxUf4u+NelEN0Mrn0XmrZ1O6IfqMT3RB/eD3vW+feanQbBBY9W+7bMzExeeuklRo4cycyZM3niiSd4/vnnWbx4MX369OH666/nueee4+677z7tdTZt2sTnn3/O4cOHSUxM5Oc//zlr165lzpw5pKenU1ZWxtChQznzzDP99C8oIgFz7ICTQEqOwcwPoXUXtyM6jkYiQaRLly6MHDkSgGuvvZbFixfTo0cP+vTpA8ANN9zAkiVLqr3ORRddRExMDLGxsXTo0IHc3FyWLl3K9OnTadq0KS1btmTKlCkB/XcRET8oLYTZV8HB7XDVG9Dx9F3G3aCRyIlqMGIIFFOLKotGjRrh9TpnmBQVFR33WkxMzA9/joyMpKysrNbXFxGXecvh3Vtg13K4/BXoEZxdujUSCSJZWVl8++23AMyePZtx48axY8cOtmzZAsDrr7/O6NGjAWdNZNWqVQC8++671V571KhRzJ8/n8LCQg4fPszChQsD9G8hInVmLXxwD2x6Dy74KwyY7nZEp6QkEkT69evHP//5TwYPHsyBAwf49a9/zSuvvMLll1/OoEGDiIiI4I477gDgwQcf5Fe/+hXnnnvucQdXncrQoUO58sorSUpK4tJLL9XZIyLBbMnjsPJlGHk3jLjd7WhOK2xawVc5HvfWzZs3H/ea2p//SD8LEZetfs3piTX4Kpj+D9c3E4JawQPaJyIiIeD7RbDwbjjjPJj6dFAkkOqETRIREQlq2Svh7RucLQFXvAaRUW5HVCNKIiIibtu/Bd64HFp0hGvegZgWbkdUY0oiFcJlbeh09DMQccHhXPjXdDARcO08aN7B7YhqRUkEaNy4MXl5eWH9S9RaS15eHo0bN3Y7FJHwUXQI3rgUju6Ha96Gdme4HVGtabMhkJCQQHZ2Nvv27XM7FFc1btyYhIQEt8MQCQ9lJfDWtZD7HVz9FsSHZhsiJREgKiqKHj16uB2GiIQLrxfS7oTtX8K056D3eLcj8pmms0RE6tunf4B178D5f4Ckq92Opk6URERE6tO3z8A3f4dht8I5v3E7mjoL6SRijOlpjHnJGDPX7VhERKq1bi4s+h30m+L0xAqBzYTVCbokYox52Riz1xiz/oTnJxljMo0xW4wx9wNYa7dZa292J1IRkVrY9iXMvwO6ng2XvBAUR9v6Q9AlEeBVYFLVJ4wxkcAzwAVAf2CGMaZ//YcmIuKDPeuco23b9YIZb0JUwymlD7okYq1dAhw44enhwJaKkUcJMAeYWu/BiYjU1sGd8K/LoHFL52jbJm3cjsivgi6JnEI8sKvK42wg3hjTzhjzDyDZGPPAqT5sjLnNGLPSGLMy3PeCiEg9qjzatqzQSSCt4t2OyO9CZZ/IyVafrLU2D7ijug9ba2cBswBSUlLCd1u6iNSfkmPw5hWQnwXXp0KHhnnEQqiMRLKBqqfTJwA5tbmAMWayMWZWQUGBXwMTEfmJ8jKYO9PpzHvpi9DtbLcjCphQSSIrgN7GmB7GmGjgKmBBbS6g80REpF5YC+//Br7/EC58DPpPcTuigAq6JGKMmQ18CyQaY7KNMTdba8uAXwCLgI3A29baDbW8rkYiIhJ4X/4VVv8Tzv0tDL/V7WgCLmyOx62UkpJiV65c6XYYItIQrXwF3rsbkq6Bqc80iM2EoONxRUQCb9MHzjRWr/Ew+W8NJoFUJ2ySiKazRCRg1s9zFtI7J8Hlr4bM0bb+EDZJRAvrIuJ3RQUw7zaYexN07A9Xvw0xzd2Oql6Fyj4REZHgsuMrpxfWoRwYfT+MuiesRiCVwiaJGGMmA5N79erldigiEsrKiuHzh+Hrp6BtD7j5Y0g46ZpzWNB0lohITe3dCC+cD1//Dc68AW5fGtYJBMJoJCIi4jOvF5b9Az59CGJawIw5kHiB21EFhWqTiDHm+lpeM8Nau9bHeEREgkuBB1J/7pyH3mcSTPk7NO/gdlRBoyYjkbFATXYkmor35QNBl0S0JiIitbb+XXjv11BeChc/CWfeGDb7P2qqJklkB05yqMlPrjKJBB1r7UJgYUpKSsPvQyAidVOYDx/cC+vehvgUuGQWtDvD7aiCUk3XRGqaepWiRSS0bV/qlO4e3g1jHoBz74FILR+fSrU/GWvtH+sjEBERV5UVw2d/hm+eVuluLdQ4vRpjOgITgCFAa5xpqzXAJ9baPYEIzp+0JiIip5S7wdl5nrsezrwJJj4M0c3cjiokVLtPxBjTzxgzF/gOuA6IAvZU/PM6YIMxZq4xpn9AI60j7RMRkZ/weuHbZ2DWWDiS65TuTn5SCaQWajISeRV4DLjGWlt84osVh0RNBV4CfubX6EREAuW40t0LKkp327sdVcipyZrIiMo/G2MirLXeE14vAd6p+BIRCX7r5jpt28vLnLbtQ29Q6a6PaltykAsoVYtIaCrMhw/ugXXvqHTXT2qURIwxg3GOpW18itezrLVd/RmYiIhfbV8C83+u0l0/q+lP8AOcEUhExRnoGTiVWRk4i/NBv1qt6iyRMHVc6W5Ple76WY26+FprE4B4oBRYCvQE/ghsBXYB/wpUgP6i6iyRMJS7AV44D775u9Oy5A513fW3Go/lrLX7jTGDrLVbK58zxhigibX2WECiExHxhdcL/34WFv8RGreCGW9B4iS3o2qQarom8jzO1FW6MWZ3ZdKw1lpACUREgkdBdkXp7hJIvBAmP6XS3QCq6UikGLgceAhoZ4zZTEVSqfhnhrV2bwDiExGpueNKd5+CoderdDfAapRErLW/rPxzRfuTjcA2IAX4OdAViAxEgCIi1apaupswDKY/r9LdelLr+jZrba4xphx4xlqbA2CMifV7ZCIiNVG1dHfsf8E5v1Hpbj3yy0/aWrvfH9cREamxsmJY/Cen91XbnnDzJ5BwpttRhZ2aLqx3s9buDHQwgaR9IiINSO4GePdW2LsBUmbChL+oaaJLarRPBNhujCkwxnxrjHkBaAKkGGOaBjA2v9I+EZEGwOt19nzMGgNH98LVb8PF/08JxEU1nc5qCyRVfCXjbDKci7NV5DtghbX2lkAEKCICOKW78++AHUsh8SKY8hQ003Ks26pNIsaYHtba7cAXFV+Vz0cDg3CSypAAxSciAmvfgfd/C94yp2V78nUq3Q0SNRmJfA50BzDGzMPpmbUGZ2/IKmBVwKITkfBWeNBJHuvfhYThcMnzziK6BI2anCfSvcrD94HBwN3AYGNMJLAWWGut/Y9ABCgiYWrbl87O88N7YOzv4Zxfq3Q3CNXqv4i19qWqj40xXXGmsjSdJSL+UVrkdN399mlo1wtu+QTiVbobrOqU1q21WUAWsNA/4YhIWNuzHubdVlG6ezNM+LMqr4Kcz0nEGLPEWjvKn8GISJjyep2Rx2d/hsat4ep3oM8Et6OSGqjLSGSk36IQkfCVv8tZ+1DpbkCkpnt4bFEmOfmFxLVuwr0TE5mWHO+362uVSkTco9LdgEpN9/DAvHUUlpYD4Mkv5IF56wD8lkhqumM9KBljmhlj/mmMecEYc43b8YhIDRUehLkzYd4t0D4Rfv6V2rYHwGOLMn9IIJUKS8t5bFGm375H0CURY8zLxpi9xpj1Jzw/yRiTaYzZYoy5v+LpS4C51tpbgSn1HqyI1E5pIfz7OXh6OHyX5pTu3vSh9n4EgLUWT37hSV/LOcXzvgjG6axXgaeB1yqfqNiP8gwwHsgGVhhjFgAJwLqKtx2fbkUkeJQWwerXYOn/wZE90P1cp/IqLtntyBocay1fbdl/2tFGXOsmfvt+dUkiARl3WmuXGGO6n/D0cGCLtXYbgDFmDjAVJ6Ek4JyuGHSjKpGwV1YM6a/D0ifgkAe6ng2XvgA9VNgZCKuzDvK/H23i39sOEN+6CTOGd2F+uoeiUu8P72kSFcm9ExP99j3rkkS+9FsU1YsHdlV5nA2MAJ4CnjbGXMRp9qoYY24DbgPo2rVrAMMUEQDKSiDjDVjyOBzKhi4jYNqz0GO01j0CYNOeQzy+6Hs+3ZhLbPNoHprcnxkjuhLTKJIRPdoFZ3WWtXas36Ko3snuOmutPQrcVN2HrbWzgFkAKSkp1s+xiUil8lJYMxu+fAwKspyjaqc8BWecp+QRADvzjvLEJ9+zYE0OzWMace/ERG48uzvNYn781T4tOd6vSeNEwbgmcjLZQJcqjxOAnNpcQIdSiQRQeRmsfQuW/C8c3AFxQ+HiJ6DXOCWPANhTUMRTn23m7RW7aBRpuGP0Gdw+qietm0bXeyw1aQX/S+B5a23xad4TA9xurX3Kn8FVsQLobYzpAXiAq4Cra3MBa+1CYGFKSsqtAYhPJDyVl8H6ufDlX+HANug8BGa8BX0mKnkEwMGjJTz35Vb++c0Oyr2WGcO7ctd5vejQsrFrMdVkJNIJ2GKM+QBnHSQTOAy0APoAY4ALqFJNVRfGmNkV14w1xmQDD1prXzLG/AJYBEQCL1trN9TyuhqJiPiLtxzWz4MvH4W8LdBxEFz1JiReqOQRAEeKy3j5q+28sGQbR0rKmJ4cz93n96FrO/cPlzXWVr9EYIyJBW7ESRaDgNbAQZw28B8Ar1lr8wIWpR+lpKTYlStXuh2GSGjyeuG7+fDFX2F/JnQYAGPuh74XQ4QKJP2tqLScN5Zl8eznW8g7WsLEAR357YRE+nRsUa9xGGNWWWtTTvZajdZErLX7gccrvkQk3Hi9sHEBfPEo7NsI7fvB5a9Cv6lKHgFQVu7l3dXZ/O3TzeQUFHFOr1jumZhIUpfWbof2E6GysF5nms4S8YHXC5vec9Y8ctdDbB+47GXoP13JIwC8XssH63fzxMffs23/UYZ0ac3jlw/h7F7B25CyRtNZP7zZOVf99ziL2p1xKqTmAA9ba4sCEqGfaTpLpAashcwP4ItHYM8653Co0ffDwEsgItLt6Bocay1ffL+PxxdlsiHnEH06NueeCYmM798REwRrTHWezqriOSARuAvYCXQDHsDZDDizLkGKSBCwFr5f5CSP3RnQpgdMfx4GXqajaQNkxY4DPPZRJst3HKBL2yb8vyuHMGVIPJER7iePmqjtXTENOMNam1/x+DtjzDJgC0GeRDSdJXIa1sKWT+Hz/4Gc1dC6G0x9FgZfqeQRIOs9BTz+cSZfZO6jfYsY/jxtIFemdCG6UWhNE9b27tgDNAXyqzzXBNjtr4ACRftERE7CWtj6mTPyyF4Brbo653oMmQGRUW5H1yBt23eEJz75nvfW7qZVkyjum9SXG8/uTpPo0JwmrG0SeR34yBjzd37cRf4fwGvGmPMq32St/cx/IYqI31kL27+Ezx+BXf+Glglw8ZOQdA00qv9dz+EgJ7+QpxZv5p1V2cQ0iuAXY3tx66ietGoS2sm6tknk9op//u6E5++o+AKwgA4HEAlWO75ypq12fg0t4uCi/3NOFGwU43ZkDVLekWKe/WIrr/97J1i47qxu/MfYXrRv0TB+3rVKItbaHoEKJNC0JiJhb+c3TvLYsRSad4ILHnNOE4xyr2VGQ3a4qJQXlm7npaXbKCwt59KhCfxqXG8S2ri/y9yfalviO8pau+Qkz8+w1s72a2QBohJfCTtZy+CL/4FtX0CzDnDub+DMGyHKfwcThbPUdM9xrdbvHtebg8dKePaLreQfK+XCQZ34zfg+9OpQv7vM/el0Jb61TSL7gJeB31trS40xrYHngWRrbR9/BBtoSiISNrJXOiOPrYuhWXsYeTekzITohvU3YTelpnt4YN66n5xjDjCqT3vunZDIoIRWLkTmX/7cJzIEeAXneNq/Aw/h9M7SGZciwcDrhS2fOOeYb/scmraD8X+CYbdAdDO3o2twHluUedIEEts8mtdmDnchovpX2zWRHGPMNGAZziFPL1lrbz/9p4KD1kSkQSsqgIw3YfkspyV7i84w7iEYdivENHc7ugbHWsvy7Qfw5Bee9PW8IyX1HJF7apVEjDFJwBvAZpwKrScrWrf/vMoGxKCkfSLSIO3f7CSOjDeh5IhzDO15v4d+U7TPIwDyjhQzb7WH2Suy2LbvKAanHPVEca3DZ72pttNZi4H7rLUvAhhjPsc553wdx588KCKB4vU66xzL/uHsMo+MhoGXwvDbIH6o29E1OF6v5dttecxensWiDXsoLbec2a0Nj1/eC6/X8uCCDcdNaTWJiuTeiYkuRly/aptEhllrt1U+qDjj/GZjzBT/hiUiP1F0qMqU1VanTHfsfzmVVs07uB1dg7P3cBFzV2Xz1opd7Mw7RqsmUVx7VjdmDO963Hke0Y0ijqvOundiYkDPNA82tV0T2WaMGQ/MANpbaycbY1KAIwGJTkRg/5aKKas3nCmrhOEw9nfOlJV2l/uV12tZumU/s5dl8enGXMq8lhE92vLrcX2YNLATjaN+2ppkWnJ8WCWNE9V2TeQu4FfAi8ClFU8X4kxpne3f0ETCmNfr9LRa9g+n2ioiypmyGnEbxJ/pdnQNzp6CIt5ZuYu3Vu4i+2AhbZtFM/OcHlw5rAtntFdhwunUdjrrbuB8a+0OY8x9Fc9twmkPH9RUnSUhofgwZMyG5c87Z5c37whjfudMWbXo6HZ0DUq51/Ll93t5c9kuPtuUi9fCyF7tuG9SXyYM6EhMo9BsiFjfaptEWgC7Kv5cWZQQBQR9PZuqsySo5W2F5S9A+r+g5DDEp8AlL0L/qZqy8jNPfiFvr9jF2yt3sbugiNjmMdw++gyuTOlC91jtpamt2iaRJcD9wMNVnvsl8LnfIhIJF14vbPsMls2CzR9DRCMYMB1G3A4JJ90cLD4qLffy2aa9zFmexRff7wPg3N7teXByf87v15GoyNA6wyOY1DaJ3AUsNMbcCrQwxmQCh4DJfo9MpKEqPgxr5sCy5yFvs9PPavR9kHITtOjkdnQNyq4Dx5izIot3Vmaz93AxHVrE8IuxvbgipQtd2qr9iz/UtjprtzFmGDAM52jcXcBya603EMGJNCgHtv04ZVV8COKGwiUvVExZNYy24MGgpMzLpxtzmb08i6Wb9xNhYExiB2YM78rYxPY00qjDr2p97qV1OjYur/gSkdOx1ulhtex55+zyiMiKKas7NGXlZ9v3H2XOiizeXZXN/iMlxLVqzN3jenNFSpew2kFe33R4skggFB+BtXOc9Y79mU4X3dH/CWfeBC07ux1dg1FcVs6iDbnMXpbFt9vyiIwwnNe3A1cP78qoPu2JjDBuh9jgKYmI+NOB7bDiRVj9OhQXQOckmP68M/rQlJXfbNl7hDnLs3h3dTYHj5WS0KYJ90zow+UpXejYUods1aewSSLaJyIBY61z4NPyWZD5oTNl1X9qxZTVMDD627A/FJWW88G63cxZvovlOw7QKMIwYUBHrhrWlXN6xRKhUYcranUo1XEfNGaJtXaUn+MJOB1KJX5TctSpslo+C/ZtgqaxToVVykxoGed2dA3Gpj2HmLN8F/NWZ3OoqIzu7Zpy1fCuXDo0ocGcUx7s/HkoVVUj6/BZkdC1Zz2kvw5rZjvneHQeAtOegwGX6LxyPzlWUsZ7a3cze3kW6Vn5REdGMHFgJ2YM78JZPdpp1BFEwmY6S6ROig7B+rnOWkfOaqf9et+LnY2BXUZoyspP1nsKmLMii7T0HA4Xl9GzfTN+f1E/LhmaQNtm2rkfjJRERE7FWsj61kkcG+ZDWSF06A8TH4HBV0Kzdm5HGJJS0z3HtU6/67xeeC3MWZHF2uwCohtFcNGgzswY3pVh3dtglKCDmpKIyIkO58KaN51NgXlbILoFDLkSkq93Dn3SLzWfpaZ7eGDeuh8OcfLkF3L/vHUAJHZswYOT+zM9OZ7WTTXqCBVKIiIA5WVOy/XVr8P3H4Eth64/g3N+AwOmQbQa8/nDox9uOu4UwEqxzWP46O5zNeoIQXVJIvqvLaEvb6sz4sh4E47scTYF/uw/IPk6aN/H7egahKPFZXy0fg+pGR72HCo66XvyjhQrgYSouiSRL/0WhUh9KjkGGxc4o46dX4GJgN4TnMTRZyJERrkdYcgrK/fy1Zb9zE/38PGGXApLy0lo04QWMY04XFz2k/erLUno8jmJWGvH+jMQkYCyFnZnOIlj3VxnN3mbHnDef0PS1drX4QfWWtZ7DjE/3cOCNTnsP1JMqyZRTB8az/TkeFK6tSEtI+e4NRGAJlGR3Dsx6M+1k1PQmog0bIUHYe07sPo1yF0HjRo7Z5MPvR66jYQIdXStq+yDx0jLyGF+uocte48QHRnBeX07MC05nrF92x93QmDlWeRVq7PunZgY1meUhzqfd6wHA2NMT+C/gFbW2stq8hntWA8DXi/sWOKMOjYuhPJiZ0Ng8nUw6HJo0trtCENeQWEpH67bzbx0D8u3HwBgWPc2TE9O4MJBnVRd1cAEZMe6cVbBbrDWvurj518GLgb2WmsHVnl+EvA3IBJ40Vr76KmuYa3dBtxsjJnrSwzSwBR4nAXy9Nchfyc0buWMOIZe5yQRqZOSMi9fZO5lfrqHxZv2UlLmpWdsM347vg/TkuN1yFOYqsuaiDXGTAFe9fESrwJPA69VPmGMiQSeAcYD2cAKY8wCnITyyAmfn2mt3evj95aGoqzEKcld/RpsXQzWCz1GOWsd/S6GKC3Y1oW1ltVZB5mf7uG9tbvJP1ZKu2bRXD28K9OT4xmc0EpVVWGurmsijYwxS3AOqPICWGv/syYftNYuMcZ0P+Hp4cCWihEGxpg5wFRr7SM4oxafGGNuA24D6Nq1q6+XkWCyL9NJHGvmwLH90CLO2dORfC207eF2dCFv+/6jzE/3kJruIevAMWIaRTBhQCcuSY7nnN6xOpNcflDXJPJ/Jzyu6wJLPM6Ru5WygRGnerMxph3wMJBsjHmgItn8hLV2FjALnDWROsYobik+4rQfWf0aZC+HiEbQZxIMvQF6ne+0YBefHThawntrc5i32kPGrnyMgbPPaMcvz+/NxAEdadFYpc/yUz4lEWPM7dba53FGByf+Ul5Sh3hONi4+5S99a20ecEeNLqzzREKTtZC9Elb/00kgJUcgtg+M/zMMuQqad3A7wpBWVFrOpxtzSU338EXmPsq8lr6dWvDABX2ZmhRPp1bqSiyn5+tI5N8V/3zPX4FUyAa6VHmcAOT448LW2oXAwpSUlFv9cT0JsKP7namq9Nedszqimjqt1odeD12Gq39VHXi9lmXbDzA/PZsP1+3hcHEZHVvGcPM5PZiWHE+/zi3dDlFCiE9JxFq7puKP64EDFYvsBmhbx3hWAL2NMT0AD3AVcHUdrymhwlsOWz93Rh2ZH4K31DkZcPJTMPASiGnhdoQh4cQuuZX7ML7PPcz8dA9p6R5yCopoFh3JpIGduWRoPGf1bKfzyMUnddonYoxZbK09v8rjT62142r42dnAGCAWyAUetNa+ZIy5EHgSpyLrZWvtwz4HePz3q5zOunXz5s3+uKT4y8GdP/avOpQNTdvB4Kuc0twO/dyOLqSc2CUXICrS0KFFDJ78IiIjDKN6xzItOZ4J/TvRJFrrSFK90+0TqWsSOe6IXGPMUmvtuT5fsB5os2GQKC2CTe8501XbKtqwnXGeM12VeCE00mY1X4x89DM8+YU/eT4q0vC7C/sxeUgcsc11pKzUTqCOxwVYa4z5G04zxtHA2jpeL2C0sB4k9qx3qqvWvgVF+dCqK4x5wOlf1bpLtR+Xkyst97Lk+30nTSAAZeWWm0aq9Fn8r0ZJxBhzD5AOpFtrD1Q+b639RcUv537ApxWL10FJC+suKsyH9e86o46c9B+Plh16HfQYo/5VPrLWsmrnQVIzPLy/djcHj5USYcB7kskFdcmVQKnpSGQScD/QxhiTDazGWQRPq/zlHKD4JFSVl8KWxbBmtrNIXl4MHQbApL/C4CugaV1rMMLX5tzDpGZ4SMvIIftgIY2jIhjfvxPTkuI4eLSE/07boC65Um9qlEQqF8uNMd2AZGAoMAr474q2JDdZa48FLEoJDdbCnrVOae66d+DoPmeRPOUm50zyuGSV5vpod0EhCzJySM3IYePuQ0RGGEb2iuU34/swYUAnmsf8+L9yo8gIdcmVelPXhfVY4E1gpbX2d36LKgBUnRVAh3Y7SWPNbNj7nTNd1WcSDJkBvcfrkCcfVXbKTc3wsGz7AayFIV1aMy0pjosHx9G+hRbIpX4ErDqr4uJ9gA+stSGxYq3qLD8pOQaZHzhluds+dxofJgxzdpEPuETTVT4qKi3n8017Sc3w8PmmfZSUO51ypybFMzUpju6xOutd6l8gq7MAsoDOfriOBDuvF7K+cUYcG9Kg5DC06uI0PhwyA2JD4u8RQafca/n3tjxS0z18tN7ZQd6+RQzXntWNaclxDIpXp1wJXjWtzsqnojoLZ1E9HdhorfUC1wBbAxWgv6jEtw7ytjrrHGvnQH4WRDeH/tOcUYdOB/SJtZYNOYdITfewcG0OuYeKaR7TiIkDOjEtOY6zz4jVDnIJCTWazjLGnAMk4SyqJwEDcBojFgIxwBXW2vcDFqUfaTqrhgoPOg0P18yBXcsAA2eMdUYcfS+CaE2r+CIr7xhpGR5SMzxs3XeUqEjDmMQOTE2KY1y/jjSO0g5yCT51ns6y1n4FfFXlgo2A/kAHYJ21NtcfgYrLfijLfbOiLLcE2veFcX90ynJbxrkdYUjKO1LM++t2k5ruYXVWPgDDe7Tl5nN66ihZCXm+NmAsI4h3p0stWAu71/xYlntsPzSNhZSbnemqzkNUluuDYyVlfLwhl9QMD0s376e8osX6fZP6MiUpjnht/pMGotokYoy5vpbXzLDWKsEEu0O7Yd3bTvKoLMtNvMCZruo1TmW5p3GqLrml5V6+2ryf1AwPH2/IpbC0nLhWjbn13J5MS46jbye1WJeGpyYjkbHU7MRCU/G+fIJwlKKFdZyy3E3vO9VVP5TlDoeLnoAB01WWWwMndsn15Bdy37trmbtqF9/tPsyBoyW0ahLF9KHxTEuKJ6VbGyK0QC4NWE2SyA6c5FCT/xMqk0jQCdveWZVluRmz4bvKstyucO5vnVFHuzPcjjCkPLYo87iWIgDFZV6+2pLHRYM7My0pntF92hPdSBVrEh5quiZS079K6a9cwSJvqzPiWPMWFFQpy02aAV3PVlmuD/YUFJ2yS64Bnrl6aP0GJBIEqk0i1to/1kcg4gfHDvxYlpu9HEwE9BwD5/+hoiy3qdsRhpyCY6V8uP7H1iOnoi65Eq78sWNd3FRWAls+cRLH9x9VlOX2g/F/gkGXqyzXB0Wl5SzeuJe0DA9fZP7YeuTu8/vQNDqSJz75Xl1yRSooiYQia2HXcudgpw3znI2BKsutk7JyL99szSMtI4dFG/ZwpLiMDi1iuO5n3ZiWFM/A+JY/tB5p3yJGXXJFKoRNEmkQ1Vn7tzhluWvfgoM7oFETZ5pq8JXObnKV5daKtZY12QWkpnt4b+1u9h8ppkVMIy4c1ImpSfGc1bPdSVuPTEuOV9IQqVDnLr6hJuTanhzdD+vnOX2rPKsAAz1Hw+CroN/FENPC7QhDztZ9R0jLyGFBhocdeceIjozgvL4dmJYcx5jEDmo9InKCQHfxFX+rbLO+9m3Y8inYcug4CCb8BQZeqnUOH+QeKmLhmhzSMnJY5ynAGDj7jHbcOaYXEwd2olUTjeJEfKEkEiy85bBjqZM4vlvg7OdoGQ9n3+X0reo4wO0IQ05BYSmL1u8hNcPDt9vysBYGJ7Ti9xf1Y/KQODq2bOx2iCIhT0nEbXvWO2sc6+bC4RyIaQkDpjrrHN3O0X6OWqo81CktI4fPMvdSUuale7um/PK83kxJiuOM9s3dDlGkQVEScUOBB9bPdUYdueshohH0Gg8TH3b6V0Vpz0FtlHst327NIy3jx0OdYpvHcM2IrkxLimdwgg51EgkUJZH6UnQINi50Fsi3LwWsc5zshY87x8k2a+d2hCHFWss6TwGp6TksXJvDvsPOoU6TBnZialIcP+vZjkaRGsWJBJqSSCCVl8LWz5yNgJkfQFkRtOkBo+9z1jnUt6rWtu8/SlqGhwUZOWzbf5ToyAjG9m3P1KR4zuuryiqR+hY2SaTe9olYC57VzjrH+ned8zmatIXka52y3IQUbQSspb2Hi1i4ZjcLMjysyXYqq87q0Y7bR/dk0oDOtGqqyioRt2ifiL8c2O4c6rT2LcjbApEx0PfCio2A50MjnV5XG4eKnMqqtIwcvtm6H6+FAXEtmZYUz8VDOtO5ldaNROqL9okESmXDw7Vv/XgOefdzYOTd0H8KNG7ldoQhpbisnM837SMtw8PiTU5lVde2TfnF2F5MSYqjVwdtrBQJNkoitVVa5DQ6XPs2bP4YvKVOw8NxDzkND1sluB1hSCn3WpZtzyMtPYcP1u/mcFEZsc2juXp4V6YmxZHUpbUqq0SCmJJITVgLO7+uaHiYBsUF0LwTjLjdma7qNEjrHKdx4nGy90zoQ++OLUhN97BwbQ65h4ppFh3JxIFOz6qRZ6iySiRUKInU1IK74HCuM001+AroMRoiVAlUnZMdJ/ubt9dggahIw+g+Hfjvi+M4v29HmkTr5ykSapREasIYuPINaNMNopu5HU1IefTDTT85TtYCrZtE8cW9Y2jdVAUHIqFMSaSmOvZ3O4KQcbiolEUbcknL8LDnUNFJ31NQWKoEItIAKImIXxSXlfNl5j7SMnL4dGMuxWVeurRtQouYRhwuLvvJ+3WcrEjDoCQiPvN6Lcu2HyAtw8MH63ZzqKiMds2iuWpYF6YkxTO0a2vSMnKOWxMBHScr0pCEdBIxxkwDLgI6AM9Yaz92N6KGz1rLhpxDpGV4WLhmN3sOFTmVVQM6MSUpjpG9YomqUllVeQKgjpMVaZhc27FujHkZuBjYa60dWOX5ScDfgEjgRWvtozW4VhvgcWvtzdW9N+RONgwSO/OOkpaRQ1qGh637jtIowjAm0elZNa6fKqtEGrJg3bH+KvA08FrlE8aYSOAZYDyQDawwxizASSiPnPD5mdbavRV//n3F58SP9h0u5v21OaRm5JCxKx+A4T3acvM5PblgYCfaNNPCuEi4cy2JWGuXGGO6n/D0cGCLtXYbgDFmDjDVWvsIzqjlOMbZyvwo8KG1dnWAQw4Lh4tK+XhDLqkZHr7e4vSs6te5JQ9c0JfJQ+K0IC4ixwm2NZF4YFeVx9nAiNO8/y5gHNDKGNPLWvuPk73JGHMbcBtA165d/RRqw/FDZdWaHD797sfKqjvHOD2r+nRUzyoROblgSyIn6x1yykUba+1TwFPVXdRaOwuYBc6aiM/RNSCVlVUL1nj4YN0eCgpLadssmiuHdWFqRWWVelaJSHWCLYlkA12qPE4Acvxx4Xo7TySIVVZWLViTw4KMHPYcKqJpRWXV1JNUVomIVCfYksgKoLcxpgfgAa4CrvbHha21C4GFKSkpt/rjeqFkZ95RFmTkkLYmhy17j/xQWfW7i/oxXpVVIlIHriURY8xsYAwQa4zJBh601r5kjPkFsAinIutla+0GP32/BjsSObFL7r0TEzmndyzvrXESR3pWPuBUVj08fSAXDuysyioR8QudbBjiTuySCxBhnO71FqeyampSHJOHxBGvyioR8UGw7hMRP/jfj37aJddroXlMI+bdebYqq0QkoMImiTSk6ayqlVU5BSfvknu0uEwJREQCLmySSKgvrFtr+W73IdIyjq+sahIV+ZORCKhLrojUj7BJIqEqK+8YaRme4yqrRvdxKqvG9evAxxty1SVXRFwTNkkklKaz9h8p5v21u0nN8PxYWdW9LX+ZNpCLBh1fWaUuuSLiJlVnBYkjxWV8vGEPqRk5fL1lP+VeS99OLZiWHK/KKhFxlaqzglRJmZcvv99HWoaHTzfmUlTqJb51E24f1ZOpSfEkdtLCuIgENyWReub1WpbvOEBaRg4frNv9Q8+qy8/swrTkOIZ2baOeVSISMsImibi5JlJZWbUgI4cFa3LYXeBUVk3o35GpyfGco55VIhKitCYSQLsOVFRWZeSwuUpl1ZSkOMb370jT6LDJ4SISwrQmUo8qK6vSMjysrqaySkQk1CmJ+EFlZVVaRg5fVamsum9SXyYP6UxCm6ZuhygiEhBKIjVwsi65Fw7qzJLv95GqyioRCWNhsyZSZWH91s2bN9f4cyfrkhsZYYiONBSWemnTNIqLB8cxNcmprIqIUGWViDQsWhPB995Zjy3K/ElvqnKvhcgIXrlxGOf0VmWViISvsEkivsrJLzzp80Wl5Yzt26GeoxERCS76K3Q1TtUNV11yRUSURKp178REmkQdfwa5uuSKiDjCZjrL1x3r6pIrInJqYVOdVSlYu/iKiASr01VnaTpLRER8piQiIiI+UxIRERGfKYmIiIjPlERERMRnYVedZYzZB+wEWgEFtfx4TT9Tk/dV957TvX6q12KB/dVG5x5ffub1de1A3g81fa8v/81P95ruh/q9RkO+H7pZa9uf9BVrbVh+AbMC9ZmavK+695zu9VO9Bqx0++fq7595fV07kPeDP+4J3Q/1f+3aXiNc74dwns5aGMDP1OR91b3ndK/7EnswCGTcdb12IO+Hmr7X1//muh8Cc+3aXiMs74ewm85qyIwxK+0pNgRJ+NH9IFUF6n4I55FIQzTL7QAkqOh+kKoCcj9oJCIiIj7TSERERHymJCIiIj5TEhEREZ8piYQJY8w0Y8wLxpg0Y8wEt+MRdxljehpjXjLGzHU7FnGHMaaZMeafFb8XrvH1OkoiIcAY87IxZq8xZv0Jz08yxmQaY7YYY+4/3TWstanW2luBG4ErAxiuBJif7odt1tqbAxup1Lda3huXAHMrfi9M8fV7KomEhleBSVWfMMZEAs8AFwD9gRnGmP7GmEHGmPdO+OpQ5aO/r/ichK5X8d/9IA3Lq9Tw3gASgF0Vbyv39RuGzfG4ocxau8QY0/2Ep4cDW6y12wCMMXOAqdbaR4CLT7yGMcYAjwIfWmtXBzhkCSB/3A/SMNXm3gCycRJJBnUYUGgkErri+fFvEeDcEKc7+P0uYBxwmTHmjkAGJq6o1f1gjGlnjPkHkGyMeSDQwYmrTnVvzAMuNcY8Rx1apWgkErrMSZ475c5Ra+1TwFOBC0dcVtv7IQ/QXybCw0nvDWvtUeCmul5cI5HQlQ10qfI4AchxKRZxn+4HOZWA3htKIqFrBdDbGNPDGBMNXAUscDkmcY/uBzmVgN4bSiIhwBgzG/gWSDTGZBtjbrbWlgG/ABYBG4G3rbUb3IxT6ofuBzkVN+4NNWAUERGfaSQiIiI+UxIRERGfKYmIiIjPlERERMRnSiIiIuIzJREREfGZkoiIiPhMSURERHymJCIiIj5TEhFxkTFmozHmiDGmpOLrSMVXP7djE6kJtT0RCQLGmJeAbdbah92ORaQ2NBIRCQ6DgfXVvkskyCiJiLjMGBOBc/a1koiEHCUREfd1xfl/cZvbgYjUlpKIiPtaAkeBaLcDEaktJRER920E1gAHjTF93Q5GpDZUnSUiIj7TSERERHymJCIiIj5TEhEREZ8piYiIiM+URERExGdKIiIi4jMlERER8ZmSiIiI+ExJREREfPb/AUlOTLIRuzuxAAAAAElFTkSuQmCC",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# 计算两个酉矩阵之间的 L1 谱距离，即 (5) 中定义的误差\n",
    "def calculate_error(U1, U2):\n",
    "    return np.abs(linalg.eig(U1 - U2)[0]).max()\n",
    "\n",
    "# 封装计算误差的函数，自由参数为每个时间块的长度 tau 和 product formula 的阶数 order\n",
    "def calculate_total_error(tau, order=1):\n",
    "    h_2 = Hamiltonian([[1, 'Z0, Z1'], [1, 'X0'], [1, 'X1']])\n",
    "    cir = Circuit()\n",
    "    # 应为总时长为 1，故 steps = int(1/tau)，注意传入的 tau 需要可以整除 1\n",
    "    construct_trotter_circuit(cir, h_2, tau=tau, steps=int(1/tau), order=order)\n",
    "    cir_U = cir.unitary_matrix().numpy()\n",
    "    U_evolve = linalg.expm(-1 * 1j * h_2.construct_h_matrix())\n",
    "    return calculate_error(cir_U, U_evolve)\n",
    "\n",
    "# 不同的参数 tau，注意它们需要可以整除 1\n",
    "taus = np.array([0.005, 0.01, 0.02, 0.05, 0.1, 0.2, 0.5, 1])\n",
    "errors = []\n",
    "\n",
    "# 记录对应不同 tau 的总误差\n",
    "for tau in taus:\n",
    "    errors.append(calculate_total_error(tau))    \n",
    "\n",
    "# 可视化结果\n",
    "plt.loglog(taus, errors, 'o-', label='error')\n",
    "plt.loglog(taus, (3 * taus**2 * (1/taus) * np.exp(3 * taus)), label='bound') # 按照 (10) 计算的一阶误差上界\n",
    "plt.legend()\n",
    "plt.xlabel(r'$\\tau$', fontsize=12)\n",
    "plt.ylabel(r'$\\Vert U_{\\rm cir} - \\exp(-iHt) \\Vert$', fontsize=12)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2fdbf168",
   "metadata": {},
   "source": [
    "下面，我们将 `order` 设置为 2，并计算二阶 product formula 的误差："
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "id": "a3db610c",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZEAAAENCAYAAADOhVhvAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAA27klEQVR4nO3dd3xUddb48c9JI4QOoabRkSo1YEcRURHEDq5tUbE87rrNXd31Wdfd9dFdffy5PpaVFcVKwIaguBYs6CqQBBAURBBIowQCoaQnc35/3ASTkEAymcmdyZz365UXmTsz957Ea8586xFVxRhjjPFGmNsBGGOMCV6WRIwxxnjNkogxxhivWRIxxhjjNUsixhhjvGZJxBhjjNci3A6gucXGxmrv3r3dDsMYY4JKenr6PlXtWvt4yCQREZkGTOvfvz9paWluh2OMMUFFRDLqOh4y3VmqulRV53To0MHtUIwxpsUImSRijDHG9yyJGGOM8VrIjIkcT1lZGdnZ2RQXF7sdiquio6OJj48nMjLS7VCMMUEiZJJI9YH12rKzs2nXrh29e/dGRJo/uACgquTl5ZGdnU2fPn3cDscYEyRCpjvreAPrxcXFdOnSJWQTCICI0KVLl5BvjRljGidkksiJhHICqWK/A2NaKE8FrHsVPB6fnzpkurOMMSYkFeXDGzfB1g+hdScYdIFPTx/USUREBgN3ArHAclV9ujmuu3htDg+/v5md+UX06tiau6YMYsaoOL9es6KigvDw8Hof10VVUVXCwqzBaUxI2rcFFsyEAztg6qM+TyAQgN1ZIvKciOSKyDe1jp8vIptFZKuI3A2gqptU9VbgSmBsc8S3eG0O97y5gZz8IhTIyS/injc3sHhtTpPO+/LLL5OcnMzIkSO55ZZbqKiooG3btvzxj39k/PjxfPXVV8c8fvTRRxk2bBjDhg3jscceA2DHjh0MHjyY22+/ndGjR5OVldX0H9oYE3y2fAj/mgRFB+C6JTDuRr9cJhBbIvOBJ4AXqw6ISDjwJDAZyAZSRWSJqm4UkenA3ZXvabL7l37Lxp2H6n1+bWY+pRU1+xWLyir47evrWbA6s873DOnVnvumDa33nJs2bWLhwoX85z//ITIykttvv51XXnmFgoIChg0bxp///GeAGo/T09N5/vnnWbVqFarK+PHjOeuss+jUqRObN2/m+eef56mnnvLiN2CMCWqq8OXj8OF90H0YzHoVOib67XIBl0RUdYWI9K51OBnYqqrbAEQkBbgY2KiqS4AlIvIu8Kq/46udQE50vCGWL19Oeno648aNA6CoqIhu3boRHh7OZZdddvR11R9/8cUXXHLJJbRp0waASy+9lM8//5zp06eTlJTEhAkTvI7HGBOkyopgyc9hwyIYMgNmPAVRbfx6yYBLIvWIA6r3y2QD40VkInAp0ApYVt+bRWQOMAcgMfH4Gfl4LQaA0x76mJz8omMD7Niahbecctz31kdVuf7663nwwQdrHH/kkUdqjHtER0cffayq9Z6vKrEYY0LIoZ2QcjXsXAtn3wtn/gaaYcZlwI2J1KOu34Sq6qeq+nNVvUVVn6zvzao6F7gfWBMVFdWkQO6aMojWkTUHtFtHhnPXlEFen3PSpEm8/vrr5ObmArB//34yMurcMPOoM888k8WLF1NYWEhBQQFvvfUWZ5xxhtcxGGOCWFYqzJ3oDKTPfBXOuqtZEggET0skG0io9jge2OlGIFWzsHw5O2vIkCH89a9/5bzzzsPj8RAZGcmTT9abEwEYPXo0N9xwA8nJyQDcdNNNjBo1ih07dngdhzEmCK19Bd75BbTvBde9Dd0GN+vl5XjdIm6pHBN5R1WHVT6OAL4HJgE5QCpwtap+29hzjx07VmvXE9m0aRODBzfvLz5Q2e/CmCBRUQ4f3AurnoY+Z8EV8yGms98uJyLpqnrMLNiA684SkQXAV8AgEckWkRtVtRy4A3gf2AQsamwCEZFpIjL34MGDvg/aGGOaU+F+eOUyJ4GMvw2uedOvCeR4Aq47S1Vn1XN8GccZPDfGmJCQuwkWzIKD2TD9CRh9ravhBFxLxF+ssqExJuh9twyePRdKC+CGd11PIBBCScS6s4wxQUsVVjzsTOHt0h/mfAqJ492OCgihJGItEWNMUCotgNd/Ch//FYZfDrP/DR38u1dfYwTcmIgxxphK+VmQMgt2fwPn3g+n3dls6z8aKmRaIoHenbVjxw6GDRvWbNf705/+xCOPPNJs1zPGNFLGl84CwgMZcPUiOP0XAZdAIISSiHVnGWOCRtrz8MJ0aN0RbloOA89zO6J6hUwSCQbl5eVcf/31jBgxgssvv5zCwkKWL1/OqFGjGD58OLNnz6akpASA3r17s2/fPgDS0tKYOHEi4LQwZs+ezcSJE+nbty+PP/740fM/8MADDBo0iHPPPZfNmzc3+89njDmBijJ499fOCvS+ZzkJpOtAt6M6rpAZExGRacC0/v37H/+F790Nuzf49uI9hsMFD53wZZs3b2bevHmcdtppzJ49m0cffZRnnnmG5cuXM3DgQK677jqefvppfvGLXxz3PN999x2ffPIJhw8fZtCgQdx2222sX7+elJQU1q5dS3l5OaNHj2bMmDE++gGNMU1WkAevXQ87PodTfw7n/gnCjl94LhCETEskGLqzEhISOO200wC45pprWL58OX369GHgQOeTyPXXX8+KFStOeJ6pU6fSqlUrYmNj6datG3v27OHzzz/nkksuISYmhvbt2zN9+nS//izGmEbY/Q38ayJkrYZL5sJ5fwmKBAIh1BJpsAa0GPxFGjFoFhERgcfj1DApLi6u8VyrVq2Ofh8eHk55eXmjz2+MaSYb34a3boXoDjD7PYgLrh6CkGmJBIPMzEy++uorABYsWMC5557Ljh072Lp1KwAvvfQSZ511FuCMiaSnpwPwxhtvnPDcZ555Jm+99RZFRUUcPnyYpUuX+umnMMY0iMcDnzwIi66D7kOdBYRBlkDAkkhAGTx4MC+88AIjRoxg//79/PKXv+T555/niiuuYPjw4YSFhXHrrbcCcN9993HnnXdyxhln1ChcVZ/Ro0dz1VVXMXLkSC677DKrPWKMm0qOwKJr4bOHYORP4Pp3oF0Pt6PySkBuBe8P1QbWb96yZUuN52z78x/Z78IYP9u/3dm+ZO93cN4DMOG2gFz/UVvQbAXvL8EwsG6MaeG2r4B/ne2Usr3mDTjl9qBIIMcTMknEGGNcowqr5sKLM6BNN7j5Y+h3jttR+YTNzjLGGH8qL4Vlv4Y1L8LAC+DSuRDd3u2ofCaok4iIzACmAt2AJ1X1A2/PpaohPwU2VMbHjGk2R3Jh4bWQtRLO+DWcfS+EtawOoID7aUTkORHJFZFvah0/X0Q2i8hWEbkbQFUXq+rNwA3AVd5eMzo6mry8vJD+I6qq5OXlER0d7XYoxrQMO9fB3LNh19dw+XMw6Y8tLoFAYLZE5gNPAC9WHRCRcOBJYDKQDaSKyBJV3Vj5knsrn/dKfHw82dnZ7N271+ugW4Lo6Gji4+PdDsOY4LfhdXj7DojpAje+Dz1Pdjsivwm4JKKqK0Skd63DycBWVd0GICIpwMUisgl4CHhPVdd4e83IyEj69Onj7duNMcbh8cDHf4EvHoXEU+DKF6FtN7ej8quASyL1iAOyqj3OBsYDPwPOBTqISH9V/WddbxaROcAcgMTERD+HaowJScWH4M2b4ft/w5gb4IKHISLK7aj8LliSSF0j3qqqjwOP1/Fc7RfOFZFdwLSoqKjg21fAGBPYtn8O7/wSDmyHCx+BcTcF/fqPhgqWJJINJFR7HA/sdCkWY4xx5GfBB/fCxsXQIRGuXQx9QmtLoWCZKpAKDBCRPiISBcwEljTmBLZi3RjjM2VF8Onf4Ilx8P37cPYf4I7VIZdAIABbIiKyAJgIxIpINnCfqs4TkTuA94Fw4DlV/baR521YUSpjjKmPKmxaAu/fCwczYeglMPkv0DHhxO9toQIuiajqrHqOLwOWNXM4xhjj2LMR/v07Z/+r7sPgkneh9+luR+W6YOnOajLrzjLGeKXoACz7LfzzdKd09tT/hTmfWQKpFHAtEX+x7ixjTKN4KmDNC7D8L1CcD2NnO2MfMZ3djiyghEwSUdWlwNKxY8fe7HYsxpgAl/EVvHeX0/JIOh0u+Bv0GOZ2VAEpZJKIMcac0MEc+PCP8M3r0D4eLn/eGTwPkTUf3giZJGLdWcaYepUVw1f/B58/CuqBs34Hp/0ComLcjizghUwSse4sY8wxVOG7d+H930N+BgyeDuf9FToluR1Z0AiZJGKMMTXs3Qzv/Q62fQJdB8N1b0PfiW5HFXRCJolYd5YxBoCifPjsb7B6LkS1gfP/BuNuhPBItyMLSiGTRKw7y5gQ5/HAupfho/uhMA/GXA/n/De0iXU7sqAWMknEGBPCMlfBe7+FXesgYQJc8wb0Gul2VC2CJRFjTMt1aBd8dB+sXwjtesGlz8Lwy23Krg9ZEjHGtDzlJbDyKfjsYfCUwRm/htN/Ba3auh1ZixMyScQG1o0JAarO1uzv3wP7t8GgqTDlr9C5r9uRtVgnTCIicl0jz7lOVdd7GY/f2MC6MS3cvi3w73tg64cQOxCueRP6T3I7qhavIS2RswFtwOuk8nX5QMAlEWNMC1V8CFb8HVY+DZExMOV/IHmOTdltJg1JIjtwkkNDRqKqkogxxviXxwNfL4CP/gQFe2HUNTDpj9C2m9uRhZSGjok0dCpDs055EJG+wB+ADqp6eXNe2xjjoux0Z5fdnHSIHwdXp0DcGLejCkknTCKqen9zBFJFRJ4DLgJyVXVYtePnA//AKY/7rKo+pKrbgBtF5PXmjNEY45LDe2D5/bDuFWjbHS55BoZfCWEhU18v4DR4dpaIdAfOA04GOuJ0W30NfKiqu30Y03zgCeDFatcOB54EJgPZQKqILFHVjT68rjEmUJWXwqp/wmd/h/JiZ4fdM38Drdq5HVnIO2H6FpHBlZ/0NwLXApHA7sp/rwW+FZHXRWSILwJS1RXA/lqHk4GtqrpNVUuBFODihp5TROaISJqIpO3du9cXYRpjmsuWD+HpU+DD/4akU+G/VsHk+y2BBIiGtETmAw8DP1HVktpPikgUzh/0ecApPo3uR3FAVrXH2cB4EekCPACMEpF7VPXBut6sqnOBuQBjx45tyEwzY4zb8n5wtmj//t/QuR9c/RoMPM/tqEwtDRkTGV/1vYiEqaqn1vOlwGuVX/5S14C9qmoecGuDTmCLDY0JDiWHYcUjzorz8CiY/GcYfxtERLkdmalDY1es7wG6+iOQE8gGEqo9jgd2uhCHMcZfVGH9Iqc87ZHdcPLVcO590K6H25GZ42hQEhGREcAmILqe5zNVNdGXgdWSCgwQkT5ADjATuNqP1zPGNKecNU6BqOzV0Gs0zHwF4se6HZVpgIbOi1sGHAGiRWSBiPxORM4XkR4i0gvo4KuARGQB8BUwSESyReRGVS0H7gDex0lmi1T128acV1WXquqcDh18FqoxpqmO7IW374B/nQMHtsPFT8JNyy2BBJEGtURUNV5EYoFM4HOcab6XAsNwWif/9FVAqjqrnuPLcJKZV2xMxJgAUlEGq/8Fnz4EZQVwyn/BWb+FaPuQF2xEteGTlUSkn6r+UO2xAK1VtdAfwfnD2LFjNS0tze0wjAlN5SWw7lX4z2NwYAf0mwTnPwRdB7odmTkBEUlX1WOaiA0dE3kGWAesFZFdVUlDnQwUFAnEWiLGuKjkCKTPh6+egMO7oNcop7b5wClWICrINaglIiKP43RdDQW6AFuoTCqV/65T1Vy/RelD1hIxphkVHXC6rVY+DUX7ofcZToGovhMteQSZJrVEVPXn1U7UHWdwexswFrgNSMTZ08oYY5w9rlY+CanzoPQIDLwAzvgVJCS7HZnxsUZXNlTVPSJSATypqjsBKgfdA5p1ZxnTDA5kwJePw5qXnLK0Qy+F038JPYad+L0mKPmkPK6q7vPFefzJKhsa40d7N8MX/89ZLChhMPJqOO1O6NLP7ciMnzV0YD1JVTP8HYwxJsjkrIEvHoVN70Bkaxh/C5xyB3SIczsy00wa2hLZLiKHcXby/QZoDYwVkY+CZXqvdWcZ4yOqkPEf+Px/4YePoVUHZ1v28bdCm4Dv2TY+1tDZWR2BkZVfoyr/HYyzMeJGIFVVb/JTjD5ls7OM8ZIqbPnASR5Zq6BNV2eR4NgbIbq929EZP/N6dpaI9FHV7cCnlV9Vx6OA4ThJ5WSfRWqMCSyeCti4GD7/f7BnA3RIgAsfcWqaR7Z2OzrjsoZ0Z30C9AYQkTdxqhl+jbM2JB1I91t0xhj3lJfC+hT44jHY/wN0GQAznobhV0B4pNvRmQDRkHoivas9fBcYAfwCGFFZtnY9sF5V/8sfARpjmllpAax5Eb78PziUAz1PhitfhJMugjBbDmZqatQUX1WdV/2xiCTidGUFfHeWDawbcwJF+ZBaubq8MA+SToPpjzv7W9nqclOPRm3A2BLYwLoxtRzZ61QRTH0WSg7BgPPg9F9Bkr+qXZtg1KRtT+o54QpVPbNpYRljXJOf5XRZrXnB2V136AxndXnPgO9YMAGkKSvWT/NZFMaY5rNvizNYvj7FeTxiJpz+C4gd4GZUJkj5ZNsTt4hIG+ApoBT4VFVfcTkkYwLXrq/h80dh49sQ0cpZ33Hqz6BjgtuRmSAWcElERJ4DLgJyVXVYtePnA//A2S34WVV9CKe64uuqulREFgKWRIypLeMrZ4Hg1g+hVXuny2rC7dC2q9uRmRYg4JIIMB94Anix6kDlVOIngclANpAqIkuAeGBD5csqmjdMYwKYKmxd7iSPzC8hpguc898w7iZo3dHt6EwzWrw2h4ff38zO/CJ6dWzNXVMGMWOU7/Y2C7gkoqorRKR3rcPJwFZV3QYgIinAxTgJJR6nMFZYfecUkTnAHIDExETfB21MoCgrgg2vwcp/Qu630D7OqSA4+jqIinE7OtPMFq/N4Z43N1BU5nzGzskv4p43nc/dvkokTUkizTlxPA7IqvY4GxgPPA48ISJTgaX1vVlV5wJzwZni68c4jXHHoZ3OFN20550Kgt2GwvQnYMRVEBHldnTGJQ+/v/loAqlSVFbBw+9vDogk8plPImiYuhKWqmoB8NMGncAWG5qWKDvNWRy4cbGzx9WgC2HCrU4ZWlsgGPJ25hc16rg3vE4iqnq2z6I4sWyg+hSSeGBnM17fmMBRUebMsFr5NOSkOYPlybdA8s3QuY/b0ZkAkF9YyotfZSDiDI/V1quj7zbODLgxkXqkAgNEpA+QA8wErnY3JGOaWUEepD/n1C0/vAs694UL/u5UEWzVzu3oTADYfbCYeV9s49VVmRSUVjCkZzu27i2gtNxz9DWtI8O5a8ogn12zIVvB/xx4RlVLjvOaVsAtqvp4UwMSkQXARCBWRLKB+1R1nojcAbyPM8X3OVX9tjHntfK4Jmjt+dZpdWx4DcqLoe/ZMO0f0H8yhNU7n8SEkG17jzB3xTbeXJNDhSrTRvTk1on9OKlH+4CYndUD2Coiy3DGQTYDh4F2wECcP/gXUG1KblOo6qx6ji8Dlnl7XhsTMUHFUwHf/9tJHjs+h4jWcPJMp3pgt8FuR2cCxDc5B3n60x9Y9s0uosLDuGpcAnPO7EtC5x9n4s0YFefTpFFbQysbxgI34CSL4UBH4ADONvDLgBdVNc9vUfqQbcBoAlrxIVj7Mqx+Bg7scKboJt8Mo6+HmM5uR2cCgKqyctt+nvp0K59v2Ue7VhFce0oSPz2tD13btfLbdevbgDFkdvGt1hK5ecuWLW6HY0xNeT/Aqmdg3StQegQSxsOE2+CkaRAeLEOXxp88HuWjTXt46tMfWJeVT2zbVtx4eh9+MiGR9tH+LxIW8kmkirVETMBQhW2fwqp/wvfvQ1gEDLvU6bKKG+12dCZAlFV4WLJuJ//87Ae25B4hoXNrbjmzH5ePiSc6svmKhPlkK/jKuur34syM6okzzTYFeEBVi30RqDEtXmkhrF/otDz2boKYWDjrtzB2NrTr4XZ0JkAUlVawMDWTf32+nZz8Ik7q0Y5/zBzJ1OE9iQgPnAkVjW0nPw0MAn4GZABJwD04K8pn+zY037KBdeO6gzlO5cD0+VB0AHoMh4ufgmGXQWS029GZAHGwsIwXv9rB81/uYH9BKeN6d+IvM4Zy9qBuSAAuIG1Ud5aI5AH9VDW/2rHOOPtaBcWon3VnmWalClmrYdXTsHEJoJWrym+HpFNtVbk5as+hYuZ9sZ1XVmZQUFrBOSd147aJ/RjXOzD+tPqqsuFuIAbIr3asNbDL+9CMaYHKS52tSFY+BTvXQqsOzkB58hzolOR2dCaA7NhXwDMrfuCN9BzKPR6mndyLW8/qx+Ce7d0OrUEam0ReAv4tIv/Hj1uR/BfwooicU/UiVf3YdyH6hnVnmWZxZC+kP++sKj+yG7r0hwsfgZNnQau2bkdnAsg3OQd5+rMfeG/DLiLCw7hyXDxzzuhHYpfg2m25sd1Z2xvwMlXVvt6H5F/WnWX8YvcGZ/v1Da9BRQn0m+S0PPpNslXl5ihVZdX2/Tz16Q+s+H4v7VpFcM0pSfz0tN50axfY42I+6c5SVdvdzZgqngrYvMxJHhlfQGQMjLoGxt8CXX23N5EJfh6Psvy7XJ7+dCtrMvOJbRvFXVMGce0pSc2yxsOfGjvF90xVXVHH8VmqusB3YRkTwAryYO1LkPYc5GdAhwSY/Gen8FPrTm5HZwJIWYWHpV87azy+33OE+E6t+cvFQ7libEKzrvHwp8aOibxRWQP9XlUtE5GOwDPAKCCgk4iNiZgmqZpllTYPvn0LKkoh8VQ47y8waKqtKg9xtTc5vHPSAIrKKpi7Yhs5+UUM6t6Ox64ayUUjAmuNhy80dkykF/A80B34P+BPOHtn/aqyQFTAszER0yglh2H9IqfVsecbiGrnbIQ4djZ0H+J2dCYA1C5BW92YpE7cPrEf55wUmGs8GsNXYyI7RWQGsAqn3Ow8Vb3FNyEaE0D2bHRaHV8vhNLDzsLAix6D4VfYLCtTQ10laAFi20bxxm2nuhBR82rsmMhI4BVgC/B74LHK+h+3VV+AaExQKi+BTUud6bmZX0J4Kxh6CYy7EeLH2cJAU0OFR1nx/V5y6ik1m3ektJkjckdjO3KXA79T1WcBROQT4HFgAzXL1xoTPA5kOGs71rwEhfugU29noHzkNdCmi9vRmQCTd6SERWnZvLo6g6z9RYQJePxcgjaQNTaJjFPVbVUPKsdBbhSR6b4Nq2FEpC/wB6CDql7uRgwmSHkqYOtHTqtjywdOK2PgBTBuNvQ9x9Z2mBpUlTWZB3h5ZSbvrt9FaYWHCX07c/f5gykuK+fexd/W6NLydQnaQNbYMZFtIjIZmAV0VdVpIjIWONLYC1fO8roIyFXVYdWOnw/8A6cM7rOq+tDx4sFJYq839vomRB3Z60zPTX8e8jOhbXc48zdO0aeO1pg2NRWUlPP2up28tDKDTbsO0a5VBFePT+Qn4xMZ0P3HuvbhYWF+LUEbyBo7JvIz4E7gWeCyysNFOF1ajR1Bmg88QbWyuiISDjwJTMbZViVVRJbgJJQHa71/tqrmNvKaJhSpQuZKSH0WNr4NnjLofYbTZXXSRRAe3Iu9jO9tzT3MyyszeSM9m8Ml5Qzu2Z7/uWQ4F4/sRZtWx/7Z9HcJ2kDW2O6sXwCTVHWHiPyu8th3ONvDN4qqrhCR3rUOJ+PsCLwNQERSgItV9UGcVotXRGQOMAcgMTHR29OYYFN8yKnbkfYc5G6EVu2dQfKxs21FuTlGWYWHD77dw0srd7By236iwsOYOqIn10xIYnRix6CfousvjU0i7YCsyu+rhpIiAV9NQ4irdn5wWiPj63uxiHQBHgBGicg9lcnmGKo6V0R2AdOioqLG+ChWE6h2f+NMz12/yCk12/NkmPY4DL8cotq4HZ0JMLsOFrFgVSYLUrPYe7iE+E6t+d35J3Hl2Hi6tPVfzfKWorFJZAVwN84f7io/Bz7xUTx1pfp6V0Oqah5wa0NOrKpLgaVjx4692cvYTCArK3a6qtLmQdYqiIiGoZfCuJucUrP2KdJU4/EoX/6Qx0srd/DRplw8qpw9qBvXTkjizIFdCQ+z+6WhGptEfgYsFZGbgXYishk4BEzzUTxV28tXiccpwdtktu1JC7V/uzNIvvZlKMyDzv3gvAdg5NUQExjFfEzgOFhYxmvpWbyyKpPt+wro3CaKOWf25erkRBI6B9cW7IGisbOzdonIOGAcTmncLGC1qnp8FE8qMEBE+gA5wEyceu7G/MhT4UzLTX0Wti4HCYNBFzjjHX0m2vRcc4z12fm89FUGS9fvpLjMw5ikTtw5aQAXDO9Bq4iWsRGiWxq1d5ZPL+ysdJ8IxAJ7gPtUdZ6IXAg8hjMj6zlVfaDek3jB9s4KYkdyYc0LkP4CHMyCtj1gzPXO9NwOoTkzxtSvuKyCpV/v5OWVGXydfZCYqHBmjIrjmvFJDOkVHFUDA0l9e2e5lkSaW7XurJu3bNnidjimoVQh4z/OosBNS53puX3Oclodgy606bnmGNv3FfDKygxeS8/mYFEZ/bu15doJSVwyOi7oa3e4KeSTSBVriQSJ4oPO5odp82DvdxDdAUb+xJmeGzvA7ehMgCmv8LD8u1xeXpnB51v2EREmTBnWg2snJDG+T2ebnusDPtnFN5jZwHqQ2PW10+rY8DqUFUCv0XDxk85Mqygb+DQ15R4uZuHqLF5dncmug8X07BDNryYPZOa4BLq1D+xysy2F1y0REVmhqmf6OB6/s5ZIACordgo9pc2D7FSIaA3DL4OxNzrTc42ppqpO+csrM/j3N7sp9yhnDIjlmglJTDqpW4sr+hQo/NESOa0J7zUG8n74cXpu0QHoMgCmPAgjZ1mZWXOMw8VlvLU2h5e+ymBL7hHaR0dww6m9+cmEJPrE2iJSt1h3lmleFeWw5X1neu4PH4OEw+CLnFZHnzNtUaA5xqZdh3hpZQaL1+ZQWFrBiPgO/P3yEUwb0YvWUTY9121N6c6qUNWg+y9o3VkuObwb1rwI6fPhUA606wVjboDR10H7nm5HZwJMSXkF//5mNy99lUFaxgFaRYQx/eReXDMhiZMTOrodXkgK+YF14wKPB7Z94nRZbX4PPOXQ92y44O8w8HwIt9vP1JS1v5BXV2eyKDWLvIJSeneJ4d6pg7l8TDwdY6LcDs/UIWT+L7burGZ0eA+se9lZFJifATFdYMJtMOan0KWf29EZly1em1Oj9savJw+kU5soXl6ZwcebcxHg3MHdufaUJE7rF0uY7WMV0JrSneVR1aCbBmHdWX5ytNUxHzYvc1odvc9wuqwGT4MI2w3VOAnknjc31KgCKDi7rMa2bcWs5ARmJSeGTGnZYOKP7qzPmvBe01LU1+oYfQPEWqvP1PTw+9/VSCDgJJBOMZF8efc5REUE3efSkOd1ElHVs30ZiAkiHg9s/xTSnq/Z6pj0R2t1mDrlHSnhjTXZ5OQX1/l8fmGZJZAgZWMipuGO5DprOta8AAd2QOvOMP5WZ6zDWh2mlqqaHQtSM/ng292UVShR4WGUVhy76bd1XwWvkEkiVpTKS1WtjvT58N27P7Y6zvlva3WYOuUeKua19GwWpmaRub+QjjGRXDuhNzOTE9i489AxYyKtI8O5a4qVKw5WXicRcXY0u15V5/suHBMw6m113GAbIJpjVHiUFVv2krI6k4825VLhUSb07cyvzxvIlKE9iI50lpQN7N4OoMbsrLumDGLGKNvKP1g1ZUxERWQ6MN934RhXeTyw/TNnXUdVqyPpdKfVcdJFEGkb2pmadh0sYlFqNovSssjJL6JLmyhuOr0PV41LoG/XtnW+Z8aoOEsaLUhTu7MiRGQFsBrwAKjqb5sclWleR/b+OMPqwHZrdZjjKq/w8MnmvSxYncmnm3PxKJwxIJbfXziYyUO62wB5iGlqEvnfWo+bvTiJiMwApgLdgCdV9YPmjiEoeTywY4Uzw+q7d51iT0mnw9l/cMY6rNVhasnaX8jC1CxeS89iz6ESurVrxW0T+3HV2EQSu9g2/aHKqyQiIreo6jPARRybOFY04jzPVZ4jV1WHVTt+PvAPnBK5z6rqQ/WdQ1UXA4tFpBPwCGBJ5HiOaXV0gvG3OCVmuw50OzoTYErLPXy0aQ8LVmfyxdZ9CHDWwK785eJEzrFt1w3et0RWVv77ThOvPx94Anix6oCIhANPApOBbCBVRJbgJJQHa71/tqrmVn5/b+X7TG1VrY70+bDpncpWx2nW6jD12r6vgJTUTN5Iz2bfkVJ6dYjmzkkDuHJsgk3HNTV4lURU9evKb78B9lcOsgvQuZHnWSEivWsdTga2quo2ABFJAS5W1QdxWi01VF73IeA9VV1T13VEZA4wByAxMbExIQa3I3th3SvODKv925xWR/IcZ6zDWh2mlqqdc1NWZ/HVtjzCw4RJJ3VjVnIiZw7sSrjtYWXq0NQxkUWqOgmOztZaCJzbxHPGAVnVHmcD44/z+p9VXrODiPRX1X/WfoGqzhWRXcC0qKioMU2ML7B5PLDjc2eGVfVWx8R7YPB0a3WYY2zNPcyC1Vm8uSabA4VlJHR2pt1ePiae7lZi1pxAU5NIZK3Hvlh5VtfHnXoH7FX1ceDxE520xS82LNjntDrS59dqdVwPXW0hl6mpuKyCd9fvIiU1k9QdB4gMF84b0oOZyQm2c65plKYmkfUi8g+czRjPAtY3PSSygYRqj+OBnU09aYvc9uTouo75P86wSjzVWh2mXpt2HSJldSZvrs3hcHE5fWLbcM8FJ3HZmHhi29ruA6bxGpREROQ3wFpgrarurzquqndU/nEeDHxU+Wm/qVKBASLSB8gBZgJX++C8Lcfh3ZVjHS9Wria3Vof5Ue16HT8/pz8IvLo6i6+z8omKCOOCYT2YOS6RCX07I1aS2DRBg+qJiMhHwEigE05LYQ3OH/u3VfVbry8usgCYCMQCe4D7VHWeiFwIPIYzI+s5VX3A22vUFrT1RDwVTk3y9PlOlUCt+LFeh60mN5XqqtdRZUC3tsxKTuSSUXF0amNVAk3j1FdPpFFFqUQkCRgFjMaZRXUWsAT4qaoW+ihWv6jWnXXzli1b3A6n4Q7mOHtYrX0JDmZBTCyM+omzrsOqBJpaTnlwObsOHrvdemzbKFL/cK61OozXfFKUSlUzgAxgceVJY4FXcdZo/L7pYfpPUA2sV5TDlg+cqblbPgD1QL9z4Ly/wqALIcI+RZofqSppGQdIWZ1VZwIByDtSagnE+EWTBtZVdZ+I3AEsI8CTSFAMrB/IcFoca1+Gw7ugbQ84/Vcw+lro1Nvt6EyAyTtSwltrc0hJzWJr7hHatoogJiqcwtJju7JsgaDxF1/UE8kEevrgPH4VsC2RijKnOmD6C86Yhwj0nwxT/xcGTIHwkCn5YhrA41H+88M+UlKzjhZ6Gp3Ykb9fPoKpw3vy4cY9Vq/DNKuGzs7Kp3J2Fs6g+lpgk6p6gJ8AP/grwBYr7wdndtW6V6EgF9rHwcS7YdQ10CHe7ehMgNl9sJjX0rJYmJZF9oGiGoWeqmp0AEe3WLd6Haa5NHR21uk4s7NGVf47FGcBYBHOAsMrVfVdv0XpAwExsF5eApuWOmMd21eAhMPA850ZVv0nQVi4O3GZgFS15frC1Ew+/s7Zcv3Ufl2YmZzIeUO6Hy30ZExz8MnsrGoniwCG4Gy/vkFV9zQ9xObhyhTfvd87iWPdq1C0HzomwujrYOQ10D7gewJNM8vMK2RR2o9brndt14orxsRz5dgEese2cTs8E6J8MjuriqqW45vV6S1XWRFsXOKs68j8EsIi4KSpztTcvmdDmG2hbX5UUl7BB9/uYWFqFl9s3UeYwMRB3fjLxQmcfVI3Im3LdROgTphEROS6Rp5znaoGXIJpttlZezY6rY6vF0DxQejcF869H0ZeDW27+ffaJuhszT1Myuos3lybw/6CUuI6tuZXkwdyxdh4enawGVUm8DWkJXI2DatYKJWvyycAWyl+nZ1VWgDfvuW0OrJTITzK2btqzPXOqnKbn2+qKSqt4J31O1mYmkVahrP54eQh3Zk5LpHT+9vmhya4NCSJ7MBJDg25s6uSSGjY9bWTODa8DiWHIHYgTPkfGDET2nRxOzoTYL7JOUhKaiZvr93J4ZJy+nZtw+8vPIlLR9vmhyZ4NXRMpKEfjVr+R6iSw07SSJ8Pu9ZBRDQMvcQZ60icYK0OU8Oh4jKWrNtJSmom3+QcolVEGFOH92RmciLjeneyVeQm6J0wiajq/c0RiL81eUwkZ41T6GnDG1BWAN2GwgUPw4grnF10jamkqqRnHCAlNYt31++iqKyCwT3b8+eLh3LxyXF0iKldhseY4BUyy6GbPCaS+qwz7jHsUhjzU4gbY62OEFV7q/WqxXz7C0p5c002C1Oz2JJ7hDZR4cwYFces5ASGx3WwVodpkbxaJxLMvF4ncmgXRLWB6Pa+D8oEjbq2Wo8KD2NIr3Zs3HmY0goPoxI7MmtcIlNH9KRNq5D5nGZaOJ+uEwlJtijQ4GwnUrtWR2mFh6+zD3LDqb2ZOS6RQT3a1fNuY1oeSyLGNFB5hYec/KK6n1S4b9rQ5g3ImAAQ1ElERAYDd+JURlyuqk+7HJJpgapvQ1If22rdhCrXkoiIPAdcBOSq6rBqx88H/oFTGvdZVX2ovnOo6ibgVhEJA/7l55BNCKnahiQlNZP/bM07ug3JRcNjeGV1JsVlnqOvta3WTShzsyUyH3gCeLHqgIiEA08Ck3FquaeKyBKchPJgrffPVtVcEZkO3F15LmOa5Ps9zjYkb63N5kBhGfGdWvPryQO5vNo2JMPjO9pW68ZUcnV2loj0Bt6paomIyCnAn1R1SuXjewBUtXYCqetc76rq1HqemwPMAUhMTByTkZHhmx/AtAgFJeW8u34XKamZrMnMJzJcOG9oD2aOS+C0frYNiTEQPLOz4oDqHc/ZwPj6XiwiE4FLcWqaLKvvdao6V0R2AdOioqLG+CRSE9RUlfXZB0lJzWLp1zs5UlJO/25tuXfqYC4ZFUcX24bEmAYJtCRS10e+eptKqvop8GlDThyw5XFNszpYWMbidU5d8k27DhEdGcZFI3oxKzmB0Ym2DYkxjRVoSSQbSKj2OB7Y6YsTN9tW8CbgqCqrtu9nYWoWyzbsoqTcw/C4Dvx1xjCmj+xF+2jbhsQYbwVaEkkFBohIHyAHmAlc7W5IJljtPVzC6+nZLErLYvu+AtpFR3Dl2ASuGpfAsLgObodnTIvg5hTfBcBEIFZEsoH7VHWeiNwBvI8zI+s5Vf3WF9ez7qzQUOFRVny/l5TUTJZvyqXcoyT37szPzunPBcN60jrK6pIb40uuJRFVnVXP8WUcZ5DcW9ad1bJlHyhkUVo2r6VlsetgMV3aRHHj6X24clwC/bq2dTs8Y1os24DRBK3Scg8fbdpDSmoWn2/ZC8AZA7oya1wCkwZ3JyrC6pIb4yvBMsXXb6wl0nJszT3CorQs3kjPJq+glF4dovn5OQO4Ymw88Z1i3A7PmJASMknExkSCW1FpBcs2OAsCU3ccICJMOHdwd65KTuDMAV0JtwWBxrgiZJKItUSCU+265H1i23D3BSdx6eg4urWLdjs8Y0JeyCQRa4kEj0PFZby9bicLq9Ulv3B4T64al8D4Pp1tQaAxASRkkogJbFV1yReszuLdDTspLvNwUo923D99KDNGWl1yYwKVJRHjqrwjJby5JoeU1Ex+2FtAm6hwLhkVb3XJjQkSIZNEbEzEPYvX5tTYOv03kwfSpV0rFqZm8cHG3ZRVKKMTO/L3y0cwdbjVJTcmmNg6EeNXi9fmcM+bG2rUJRecXTU7xURy6eh4rhqXwMDuVpfcmEAW8utEjDv+/v53NRII/JhAVv5+Eq0ibBsSY4KZJRHjFxl5BSxMzWJnfnGdz+cXllkCMaYFCJkkYmMi/ldXXfLoiDCKyz3HvLZXx9YuRGiM8bWQSSK2TsR/tuY6dcnfWOPUJY/r2JpfTR7IFWPjWbVt/zFjIq0jw7lryiAXIzbG+ErIJBHjW3VtQzJ5SHdmJidyev/Yo9uQzBgVB1BjdtZdUwYdPW6MCW6WREyjfLvzICmrs1i8LofDxT9uQ3LZ6Hi6tqu7LvmMUXGWNIxpoSyJmBM6UlLOknU7SUnNZH32QaIiwrhwWA9mJifaNiTGhLigTyIi0gZYgVMZ8R2342kpVJV1WfmkrM5i6fqdFJZWMKh7O+6bNoRLRsXRMSbK7RCNMQHAzfK4zwEXAbmqOqza8fOBf+CUx31WVR86wal+ByzyW6Ah5mBhGW+tzSYlNYvvdh+mdWQ4007uyczkREYldLRWhzGmBjdbIvOBJ4AXqw6ISDjwJDAZyAZSRWQJTkJ5sNb7ZwMjgI2A7QneBKrK6u37SUnNYtmGXZSUexge14EHLhnG9JN70S7aNj80xtTNzRrrK0Skd63DycBWVd0GICIpwMWq+iBOq6UGETkbaAMMAYpEZJmqHrMoQUTmAHMAEhMTffpzBLO8IyW8scZpdWzbW0C7VhFcMTaemeMSGRbXwe3wjDFBINDGROKArGqPs4Hx9b1YVf8AICI3APvqSiCVr5srIruAaVFRUWN8F27w8XiUL7buIyU1kw837qGsQhmT1ImHL+/H1BE9iYkKtFvCGBPIAu0vRl0d7ifcIVJV5zfgNSG92HD3wWJeS8tiYVoW2QeK6BQTyXWn9LbND40xTRJoSSQbSKj2OB7Y6YsTh+K2J+UVHj7dvJeU1Ew+/i4Xj8Kp/brw2/NPYsrQ7rZ3lTGmyQItiaQCA0SkD5ADzASudjekwFe7Xsfs03uTX1jGa2nZ7D5UTGzbVtxyVj+uGptA79g2bodrjGlBXKsnIiILgIlALLAHZ53HPBG5EHgMZ0bWc6r6gC+v29LqidRVr6PKxEFdmTkukUmDuxEZHuZCdMaYliLg6omo6qx6ji8Dlvn6ei2hO+tISTmZeYVk7i8gI6+QjP2FvJGeTUkdu+R2b9+K+T9NdiFKY0woCbTuLL8JhoF1VWXv4RIy9xceTRKZeQWV/xaSV1Ba4/WdYiLrTCAAuYdKmiNkY0yIC5kkEigtkdJyDzn5RWTkFZBZmRyqkkTm/sIa3VJhAj07tCaxcwyTh3QnsUsMSZ3bkNQlhoTOMXRoHclpD31MTn7RMdexeh3GmOYQMkmkKS2R2gPXJ9rK/HBxGRmVSSGjsvup6vud+UV4qg1DRUeGkdg5hsTObTh9QKzzfZcYkjrHEN8phqiI449l3DVlkNXrMMa4JmSSiLctkdoD1zn5Rdzz5nryi0oZ3KO905qo1f10oLCsxjk6t4kisXMMY5I6cemoOBK7tCGxcwxJXWLo1q5Vk/ajsnodxhg3uTY7yy2NnZ1VX3dRdWHidB8ldXFaFEmVLYmEykRhe08ZY4JdwM3OChY7j5NAXpidTFLnGHp1bH3CbidjjGmJLImcQK+OretsicR1bM1ZA7u6EJExxgSOkPn4LCLTRGTuwYMHG/W+u6YMonVkze1BbODaGGMcIZNEVHWpqs7p0KFxW5zPGBXHg5cOJ65jawSnBfLgpcNt4NoYY7DurAaZMSrOkoYxxtQhZFoixhhjfC9kkoi3YyLGGGPqFzJJxNsxEWOMMfULmSRijDHG9yyJGGOM8VrIbXsiInuBDKAD0NgBksa850Svbcrz9T0XC+xrUHTu8OZ33pznbuw5fHk/nOg1dj80/7n9+TeiqffD8Z731/2QpKrHrrBW1ZD8Aub68z0nem1Tnq/vOSDN7d+rr3/nzXnuxp7Dl/eDt//N7X4InPuhMe9p6v1wgv/uzXo/hHJ31lI/v+dEr23K897EHgj8Gbcvzt3Yc/jyfjjRa+x+aP5z+/NvRFPvh+M936z3Q8h1Z7VkIpKmdeyyaUKT3Q+mOn/dD6HcEmmJ5rodgAkodj+Y6vxyP1hLxBhjjNesJWKMMcZrlkSMMcZ4zZKIMcYYr1kSCREiMkNE/iUib4vIeW7HY9wlIn1FZJ6IvO52LMYdItJGRF6o/LvwE2/PY0kkCIjIcyKSKyLf1Dp+vohsFpGtInL38c6hqotV9WbgBuAqP4Zr/MxH98M2Vb3Rv5Ga5tbIe+NS4PXKvwvTvb2mJZHgMB84v/oBEQkHngQuAIYAs0RkiIgMF5F3an11q/bWeyvfZ4LXfHx3P5iWZT4NvDeAeCCr8mUV3l7QKhsGAVVdISK9ax1OBraq6jYAEUkBLlbVB4GLap9DRAR4CHhPVdf4OWTjR764H0zL1Jh7A8jGSSTraEKDwloiwSuOHz9FgHNDHK+G78+Ac4HLReRWfwZmXNGo+0FEuojIP4FRInKPv4Mzrqrv3ngTuExEnqYJW6VYSyR4SR3H6l05qqqPA4/7LxzjssbeD3mAfZgIDXXeG6paAPy0qSe3lkjwygYSqj2OB3a6FItxn90Ppj5+vTcsiQSvVGCAiPQRkShgJrDE5ZiMe+x+MPXx671hSSQIiMgC4CtgkIhki8iNqloO3AG8D2wCFqnqt27GaZqH3Q+mPm7cG7YBozHGGK9ZS8QYY4zXLIkYY4zxmiURY4wxXrMkYowxxmuWRIwxxnjNkogxxhivWRIxxhjjNUsixhhjvGZJxBhjjNcsiRjjMhHZJCJHRKS08utI5ddgt2Mz5kRs2xNjAoSIzAO2qeoDbsdiTENZS8SYwDEC+OaErzImgFgSMSYAiEgYTv1rSyImqFgSMSYwJOL8/7jN7UCMaQxLIsYEhvZAARDldiDGNIYlEWMCwybga+CAiJzkdjDGNJTNzjLGGOM1a4kYY4zxmiURY4wxXrMkYowxxmuWRIwxxnjNkogxxhivWRIxxhjjNUsixhhjvGZJxBhjjNcsiRhjjPHa/wfo+P3dGGSiYgAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "taus = np.array([0.005, 0.01, 0.02, 0.05, 0.1, 0.2, 0.5, 1])\n",
    "errors = []\n",
    "\n",
    "for tau in taus:\n",
    "    errors.append(calculate_total_error(tau, order=2)) # 通过 order=2 指定 Suzuki 分解的阶数    \n",
    "\n",
    "plt.loglog(taus, errors, 'o-', label='error')\n",
    "plt.loglog(taus, (2 * taus * 3 )**3 / 3 * (1/taus) * np.exp(3 * taus), '-', label='bound') # 按照 (12) 计算的二阶误差上界\n",
    "plt.legend()\n",
    "plt.xlabel(r'$\\tau$', fontsize=12)\n",
    "plt.ylabel(r'$\\Vert U_{\\rm cir} - \\exp(-iHt) \\Vert$', fontsize=12)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8c468fd4",
   "metadata": {},
   "source": [
    "可以看到，实际计算出的模拟误差都小于其理论上界，说明这样的结果是符合预期的。\n",
    "\n",
    "## 小结\n",
    "\n",
    "本教程介绍了 Paddle Quantum 中构建模拟时间演化电路的功能，并对其背后的理论基础——Suzuki product formula 进行了介绍。利用该功能，用户可以通过搭建对于自定义哈密顿量的任意阶 product formula 电路来模拟不同物理系统的时间演化过程。\n",
    "\n",
    "量子模拟本身是一个比较宽泛的话题，其应用也十分广泛：凝聚态物理中的多体局域化、时间晶体、高温超导、拓扑序的研究；量子化学中的分子动力学模拟、反应模拟；高能物理中的场论模拟；乃至核物理和宇宙学中的相关研究。本教程中介绍的 Suzuki product formula 与基于通用量子计算机的数字量子模拟只是量子模拟的一部分，许多不基于通用量子计算机的量子模拟器，例如在冷原子、超导、离子阱以及光子等平台上做的模拟量子模拟（analogue quantum simulation）也是一个十分重要的方向。鉴于篇幅，本教程中无法进一步对这些内容逐一展开介绍，感兴趣的读者可以阅读 14 年的这篇综述文章 [8]。对于一些更新的结果，也可以参考这篇中文综述 [9]。\n",
    "\n",
    "在后续的教程 [模拟一维海森堡链的自旋动力学](./SimulateHeisenberg_CN.ipynb) 中，我们以凝聚态物理中的自旋模型为例，进一步地展示了如何利用本文中介绍的内容来进行量子多体模型的动力学模拟。同时，该教程也将介绍如何搭建不基于 Suzuki 分解的自定义时间演化电路。"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "01705a8d",
   "metadata": {},
   "source": [
    "---\n",
    "\n",
    "## 参考资料\n",
    "\n",
    "[1] Feynman, R. P. \"Simulating physics with computers.\" International Journal of Theoretical Physics 21.6 (1982).\n",
    " \n",
    "[2] Lloyd, Seth. \"Universal quantum simulators.\" [Science (1996): 1073-1078](https://www.jstor.org/stable/2899535).\n",
    "\n",
    "[3] Childs, Andrew M., et al. \"Toward the first quantum simulation with quantum speedup.\" [Proceedings of the National Academy of Sciences 115.38 (2018): 9456-9461](https://www.pnas.org/content/115/38/9456.short).\n",
    "\n",
    "[4] Nielsen, Michael A., and Isaac Chuang. \"Quantum computation and quantum information.\" (2002): 558-559.\n",
    "\n",
    "[5] Tran, Minh C., et al. \"Destructive error interference in product-formula lattice simulation.\" [Physical Review Letters 124.22 (2020): 220502](https://journals.aps.org/prl/abstract/10.1103/PhysRevLett.124.220502).\n",
    "\n",
    "[6] Childs, Andrew M., and Yuan Su. \"Nearly optimal lattice simulation by product formulas.\" [Physical Review Letters 123.5 (2019): 050503](https://journals.aps.org/prl/abstract/10.1103/PhysRevLett.123.050503).\n",
    "\n",
    "[7] Campbell, Earl. \"Random compiler for fast Hamiltonian simulation.\" [Physical Review Letters 123.7 (2019): 070503](https://journals.aps.org/prl/abstract/10.1103/PhysRevLett.123.070503).\n",
    "\n",
    "[8] Georgescu, Iulia M., Sahel Ashhab, and Franco Nori. \"Quantum simulation.\" [Reviews of Modern Physics 86.1 (2014): 153](https://journals.aps.org/rmp/abstract/10.1103/RevModPhys.86.153).\n",
    "\n",
    "[9] 范桁. \"量子计算与量子模拟.\" [物理学报 67.12(2018):16-25](http://wulixb.iphy.ac.cn/article/id/72211)."
   ]
  }
 ],
 "metadata": {
  "interpreter": {
   "hash": "3b61f83e8397e1c9fcea57a3d9915794102e67724879b24295f8014f41a14d85"
  },
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.13"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}
